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OCEAN  AND  SHIP  WAVE  MODIFICATION 
BY  A  SURFACE  WAKE  FLOW  PATTERN 


INTRODUCTION 

Ocean  surface  waves  are  altered  by  interactions  with  mean  and  fluctuating  surface  currents. 
These  waves  may  be  ambient  or  ship-generated,  or  may  be  a  combined  system  produced  by  both 
sources  of  generation.  An  important  wave-current  interaction  is  caused  by  the  interaction  of  the 
waves  with  the  surface  wake  produced  by  the  passage  of  a  surface  ship.  Remote  observations  of  a 
ship  wake  using  synthetic  aperture  radar  (SAR)  show  distinct  features  such  as  a  dark  trailing  center- 
line  region  and  bright-line  images  aligned  at  some  angle  to  the  ship’s  path.  The  dark  region  of  rela¬ 
tively  low  radar  backscatter  is  usually  associated  with  a  “dead  water”  region  which  is  relatively 
devoid  of  short  wave  components,  while  the  bright  line  feature  suggests  a  region  of  enhanced  radar 
return  within  the  apparent  angular  confines  of  the  ship’s  Kelvin  wave  pattern. 

In  this  report  the  interactions  of  ship-generated  Kelvin  waves  and  ambient  waves  with  the 
momentum  wake  of  a  surface  vessel  are  studied  to  assess  the  effectiveness  of  the  momentum  wake  in 
the  generation  of  these  remotely-observed  features.  Other  factors  such  as  the  limiting  wave  steepness 
at  the  onset  of  wave  breaking  also  are  considered. 

A  subsequent  section  of  this  report  will  briefly  survey  related  hydrodynamic  studies  which  are 
available  in  the  literature.  A  fundamental  point  from  all  of  these  studies  is  the  importance  of  the 
wave  action,  or  the  energy  divided  by  the  wave  intrinsic  frequency,  as  the  appropriate  conserved 
quantity  for  wave  motions  in  the  presence  of  surface  currents.  It  is  the  wave  action,  and  not  the  wave 
energy,  which  must  be  considered  when  surface  gravity  waves  are  modified  by  their  interactions  with 
currents.  The  basic  theory  to  be  described  will  show  that  the  interaction  between  the  waves  and  the 
ship  wake  is  described  by  characteristic  equations  for  the  kinematics  and  the  conservation  of  wave 
action  spectral  density  for  continuous  spectra,  defined  as  the  energy  per  unit  frequency.  Solution  of 
these  equations  gives  the  changes  in  the  paths  of  the  waves  propagating  through  the  ship  wake,  the 
changes  in  their  wavenumbers,  and  the  changes  in  the  amplitudes  of  the  spectral  components. 

The  results  of  application  of  the  theory  include  predictions  of  Kelvin  waves  for  a  high-speed 
destroyer,  ambient  wave  modifications,  ship  wave  modifications,  and  a  brief  comparison  with  experi¬ 
ments.  The  modification  of  an  ambient  wave  field  by  a  surface  current  pattern  is  demonstrated  by  the 
spatial  variation  of  the  spectral  composition  and  the  surface  displacement  of  the  wave  pattern.  The 
effects  of  wave  reflection,  refraction  and  breaking  are  accounted  for  in  the  model.  It  is  shown  that 
the  near- field  centerline  region  of  the  ship  wake  is  remarkably  free  of  short  wave  components,  but  the 
downstream  extent  of  this  region  is  not  sufficient  to  explain  the  long  distances  over  which  the  dark 
stem  region  of  radar  return  persists.  The  turbulent  wake  of  a  ship  alters  the  Kelvin  wave  pattern  gen¬ 
erated  by  the  ship’s  passage.  This  interaction  leads  directly  to  a  narrower  than  usual  wedge-shaped 
region  within  which  the  divergent  waves  are  dominant.  The  angular  directions  of  these  waves  are 
shifted  so  that  they  travel  at  angles  more  nearly  normal  to  the  ship  track.  The  lengths  of  the  waves 
also  are  reduced,  which  may  cause  them  to  steepen. 
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SAR  IMAGERY  OF  SHIP  WAKES 


One  of  the  striking  features  of  the  remote  observation  of  surface  ship  wakes  using  synthetic 
aperture  radar  (SAR)  imagery  is  the  apparent  variability  of  the  wake  angles.  The  usual  Kelvin  wave 
pattern  generated  by  the  steady  motion  of  a  ship  is  included  within  an  angular  region  between  ±19.5 
degrees  from  the  ship  track.  However,  many  of  the  SEASAT  satellite-based  SAR  images  produced 
wake  half-angles  in  the  ranges  from  14  to  17  degrees  and  as  small  as  1  to  7  degrees.  A  typical  exam¬ 
ple  of  a  SAR  image  of  a  surface  ship  wake  is  shown  in  Fig.  1. 

A  program  of  sea  trials  was  conducted  jointly  by  the  Naval  Ocean  Systems  Center  (NOSC)  and 
the  Canadian  Defense  Research  Establishment  Pacific  (DREP)  to  determine  the  detectability  of  wakes 
by  SAR  ocean  wave  pattern  imaging  (Hammond  et  al.,  1985).  Ship-generated  waves  of  known 
wavelength  and  direction  were  observed  with  an  airborne  L-band  SAR  to  study,  under  selected  sea 
conditions,  the  hydrodynamic  contributions  to  the  detectability  of  the  wake.  The  results  of  the  study 
led  to  the  initial  conclusion  that  the  narrow-v  angles  of  the  ship  wave  pattern  were  produced  by  an 
interaction  between  the  ship-generated  waves  and  longer  waves  present  in  the  test  area.  Lyden  et  al. 
(1985)  analyzed  much  of  the  narrow-v  wake  data  from  these  experiments  and  concluded  that  the 
bright  images  were  produced  by  ship-generated  Bragg  waves  moving  away  from  the  ship  track. 
There  were  only  light  winds  at  the  field  test  site  during  these  experiments. 

A  review  of  the  available  synthetic  aperture  radar  images  of  the  narrow-v  wakes  available 

through  mid-1985  was  made  by  Witting  and  Vaglio-Laurin  (1985).  Among  the  images  discussed 

were  those  taken  from  SEASAT  and  the  airborne  SAR  images  obtained  during  the  DREP  experi¬ 
ments.  A  proposed  mechanism  for  the  generation  of  the  narrow-v  wakes  was  based  upon  the 

development  of  turbulent  eddies  with  length  scales  on  the  order  of  the  Bragg  wavelength  in  the  boun¬ 
dary  layer  and  near  wake  regions  of  the  ship.  Wave  breaking  also  was  thought  to  play  an  important 
role  in  enhancing  the  radar  image  of  the  wake.  No  actual  predictions  of  the  observed  wake  confi¬ 
gurations  were  made  as  part  of  the  study. 

It  has  been  concluded  by  Munk  et  al.  (1987)  that  the  SAR  images  of  the  narrow-v  wakes  are 
caused  by  intermittency  or  unsteadiness  in  the  divergent  waves  of  the  Kelvin  wake  producing 
wavelengths  in  the  Bragg  scattering  range  of  the  radar.  This  generally  supports  the  previous 
hypothesis  of  Lyden  et  al.  (1985).  High  resolution  photographs  of  ship  wakes  taken  during  a  1985 
Shuttle  mission  also  were  analyzed  by  Munk  et  al.  In  contrast  to  the  radar  images  it  was  concluded 
that  the  optical  images  were  due  to  interference  fringes  resulting  from  a  complex  geometry  of  sources 
in  the  ship  wake.  Low  wind  conditions  prevailed  in  the  vicinity  of  the  wake  detection  and  observa¬ 
tion. 


The  most  recent  remote  observations  of  v  wake  patterns  were  made  during  a  series  of  experi¬ 
ments  conducted  with  the  Shuttle  Imaging  Radar-B  (SIR-B)  during  a  mission  of  the  space  shuttle 
Challenger  in  October  1984  (A.M.  Reed,  1986).  The  SAR  was  an  L-band  instrument  which  was 
deployed  during  several  orbits  of  the  shuttle  to  record  the  radar  images  of  surface  ship  wakes.  No 
narrow-v  wakes  with  half-angles  of  1  to  7  degrees  were  observed.  Rather,  the  half-angles  of  the  v 
wake  patterns  were  clustered,  with  some  variability,  in  the  range  of  14  to  16  degrees.  All  of  the  v 
wakes  were  observed  when  the  ships  were  moving  in  the  azimuthal  direction  or  at  45  degrees  to 
azimuth  relative  to  the  shuttle  orbit  direction.  In  all  of  these  experiments  the  azimuthal  direction  of 
ship  movement  was  parallel  to  the  shuttle  direction. 

In  contrast  to  the  bright  v  wake  images  which  are  often  observed  using  synthetic  aperture  radar, 
there  also  is  a  relatively  narrow,  stem-like  dark  image  along  the  ship  track  which  commonly  persists 
for  many  kilometers  downstream.  The  dark  portion  of  the  radar  image  suggests  that  this  region  of  the 
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Fig.  1  —  A  digital  ground  range  SAR  image  of  the  narrow-v  wake  of  the  USS  QUAPAW  during  the 
1983  Georgia  Strait  experiments;  from  Lyden  et  al.  (1985). 
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surface  is  relatively  smooth  and  free  of  short- wavelength  scatterers.  This  region  has  been  described 
as  the  turbulent  wake  or  scar  of  the  ship,  but  the  recent  computations  by  Swean  (1987)  to  simulate  the 
wake  of  a  twin-screw  high  speed  destroyer  from  model  test  data  suggest  strongly  that  the  surface 
features  of  the  classical  turbulent  ship  wake  do  not  persist  sufficiently  far  downstream  to  explain  the 
SAR  images.  These  computations  are  summarized  in  a  later  section  of  this  report.  Another  sugges¬ 
tion  based  upon  the  results  of  a  recent  NRL  field  experiment  is  that  the  persistence  of  the  dark  scar 
region  is  due  to  the  redistribution  of  surfactant  materials  at  and  near  the  water  surface  as  a  result  of 
the  ship’s  passage  (Kaiser  et  al.,  1987). 

RELATED  HYDRODYNAMIC  STUDIES 

Among  the  earliest  and  most  comprehensive  studies  of  the  interaction  between  linear  surface 
waves  and  currents  is  the  work  of  Longuet-Higgins  and  Stewart  (1961,  1964).  Because  of  the 
assumptions  inherent  in  these  studies,  they  were  limited  to  the  propagation  of  short  waves.  Thus  the 
characteristic  length  of  the  wave  train  is  much  less  than  the  length  scale  of  the  flow  field  over  which 
it  progresses.  Peregrine  (1971)  studied  the  interaction  between  a  ship’s  wake  and  its  wave  system 
using  a  formulation  based  upon  the  work  of  Longuet-Higgins  and  Stewart.  Except  for  very  weak  sur¬ 
face  wakes  the  included  half-angle  of  the  wave  pattern  was  found  to  be  less  than  that  of  the  well- 
known  Kelvin  ship  wave  pattern  (Ursell,  1960).  Under  certain  limiting  conditions  the  radiated  steady 
ship  waves  at  infinity  depend  only  on  the  maximum  value  of  the  surface  current  in  the  wake. 

More  recently,  Peregrine  and  Thomas  (1979)  have  considered  the  same  problem  for  finite- 
amplitude,  symmetric  waves.  This  approach  is  based  upon  the  conservation  of  wave  action  as  dis¬ 
cussed  by  Lighthill  (1979)  and  Crapper  (1984),  and  in  the  related  work  of  James  et  al.  (1977), 
Cooper  and  Griffin  (1986),  and  Cooper,  Griffin  and  Wang  (1987).  The  finite-amplitude  wave  proper¬ 
ties  were  computed  by  Peregrine  and  Thomas  using  the  method  of  Longuet-Higgins  (1975)  for 
predicting  the  energy  and  wave  action  densities  and  fluxes.  These  studies  show  that  it  is  the  wave 
action  and  not  the  wave  energy  which  should  be  considered  when  the  gravity  wave  field  is  modified 
by  its  interaction  with  a  surface  flow  field. 

Savitsky  (1970)  conducted  an  experimental  study  of  the  interaction  between  deep  water  gravity 
waves  and  a  near-surface  turbulent  flow  field.  The  flow  consisted  of  mean  and  random  turbulent 
wake  components  which  were  generated  by  the  steady  motion  of  a  grid  in  the  same  direction  as  the 
wave  motion.  The  presence  of  the  mean  velocity  gradients  in  and  near  the  surface  plane  resulted  in 
wave  refraction,  diffraction  and  interference  effects  which  reduced  the  wave  heights  in  the  wake 
region  to  approximately  10  percent  of  their  original  value.  Outside  of  the  wake  the  resultant  local 
wave  heights  were  increased  by  as  much  as  75  percent  above  their  original  values. 

Studies  of  irregular  or  random  wave  field  interactions  as  considered  here  were  conducted  by 
numerous  workers,  e.g.  Karlsson  (1969),  Tayfun  et  al.  (1976),  James  et  al.  (1977),  Sakai  et  al. 
(1983),  and  Hedges  et  al.  (1985).  The  results  of  Tayfun  et  al.  and  James  et  al.  are  most  relevant  here 
as  they  form  the  basis  for  the  work  to  be  outlined  in  this  report.  For  more  complete  details  of  the 
model,  the  reader  is  referred  to  Cooper  (1986).  The  mathematical  foundations  for  the  overall 
approach  taken  here  are  discussed  by  Skop  (1984). 

MODEL  TEST  DATA 

The  model  test  data  employed  in  this  study  were  obtained  from  towing  basin  tests  of  a  twin- 
screw  high  speed  destroyer  model  (Lindenmuth,  1986).  The  model  and  test  conditions  as  they  relate 
to  the  present  study  are  summarized  in  Table  1.  Mean  velocities,  rms  turbulence  velocities,  and  the 
Reynolds  stresses  were  obtained  at  several  planes  aft  of  the  model  using  a  laser-Doppler  anemometer 
(LDA)  and  a  hot-film  anemometer  array. 
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Table  1  —  Twin-Screw  High  Speed 
Destroyer  Ship  Parameters  (1/25  scale) 


I 


Model 

High  Speed  Destroyer 

DTNSRDC  No. 

5415-01 

Length* 

18.8  ^ 

Beam 

29.9  in 

Draft 

9.8  in 

Prop  Dia. 

8.16  in 

Block  Coef. 

0.506 

Tow  Speed 

4.0  knots 

Froude  No. 

0.28 

*Length  between  forward  and  aft  perpendiculars. 


These  data  subsequently  were  employed  in  a  numerical  simulation  of  the  downstream  evolution 
of  the  model’s  turbulent  wake  (Swean,  1987).  The  initial  plane  of  the  calculation  was  located  four 
ship  beams  downstream  of  the  stem  of  the  model.  A  contour  plot  in  the  y-z  plane  of  the  downstream 
(x-direction)  mean  velocity  at  this  location  is  shown  in  Fig.  2,  for  the  experiments  with  outboard  pro- 
pellor  rotation.  These  data  represent  an  interpolation  of  the  experimental  data  directly  to  the  compu¬ 
tational  grid.  The  calculation  extended  downstream  to  approximately  twelve  ship  beams  aft.  The 
velocities  in  the  initial  plane  at  three  typical  depths  below  the  free  surface  are  plotted  in  Fig.  3.  The 
velocity  here  is  measured  relative  to  still  water,  so  that 


Am  _  U(y,z) 
U  U 


(1) 


where  U  is  the  ship  velocity  (tow  speed).  The  negative  values  of  Am  are  in  the  direction  of  the  ship 
motion  and  are  representative  of  a  drag  wake.  As  the  depth  increases,  there  is  a  change  of  sign  and 
Am  becomes  positive  due  to  the  propellor  thrust. 

A  contour  plot  of  the  surface  (x-y)  plane  distribution  of  mean  velocity  Am  in  the  x  direction  is 
shown  in  Fig.  4a,  where  the  downstream  and  cross-stream  distances  are  scaled  by  the  ship  beam  B. 
The  computed  results  show  that  the  mean  surface  features  of  the  wake  are  greatly  diminished  at  x  /  B 
=  12,  the  downstream  limit  of  the  experimental  domain.  This  distance  corresponds  to  approximately 
1.6  ship  lengths,  so  it  is  clear  that  the  free  surface  mean  velocity  field  in  the  wake  of  a  ship  has 
greatly  diminished  after  but  a  few  ship  lengths  aft  of  the  vessel.  The  calculation  recently  has  been 
extended  to  1000  ship  beams  downstream,  and  the  results  are  shown  in  Fig.  4b.  In  the  far  wake  the 
propellor  thrust  has  broached  the  surface  and  has  displaced  the  drag  wake  outboard. 

An  estimate  of  the  width  of  the  early  turbulent  wake  at  the  surface  can  be  made  from  the  results 
of  the  calculation.  If  the  contour  of 
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Fig.  2  —  The  initial  (y  -  z)  mean  velocity  plane  for  the  calculation  of  the  turbulent  wake  of  the 
Model  5415  high  speed  destroyer  model. 
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Fig.  3  —  Measured  mean  velocities  Au/U  at  three  typical  depths  in  the  initial  (y  -  z)  plane  of  the 
Model  5415  wake.  (1  ft  =  0.305  m) 
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Fig.  4a  —  Contour  plots  of  the  surface  (x,y)  plane  distribution  of  mean  velocity  A u/U  from  the  cal¬ 
culation  of  the  Model  5415  wake;  from  Swean  (1987). 

The  early  wake,  0  <  (*  —  x0)/B  <  12. 
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Fig.  4b  —  Contour  plots  of  the  surface  (x,y)  plane  distribution  of  mean  velocity  A u/U  from  the  cal¬ 
culation  of  the  Model  5415  wake. 

The  late  wake,  0  <  (x  —  x0)/B  <  1000. 
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which  is  shown  in  Fig.  4a,  is  chosen  arbitrarily  as  the  cross-stream  limit  of  the  wake,  then  the  half¬ 
width  d  /  B  —0.75.  This  is  approximately  half  of  the  wake  width  which  typically  is  measured  from 
aerial  photographs  of  the  visible  “white  water”  wake  of  a  surface  ship  (Peltzer,  1984).  A  compar¬ 
able  representation  of  the  wake  half-width  d  at  500  ship  beams  (500  B)  downstream  in  Fig.  4b  is  d/B 
=  1.25  when  the  Au  =  0  contour  is  employed  there  to  define  the  cross-track  extent  of  the  surface 
wake.  This  is  equivalent  to  66  ship  lengths  for  the  destroyer  model  with  B/L  =  0.133.  A  typical 
example  of  the  visible  white  water  wake  from  the  study  by  Peltzer  is  shown  in  Fig.  5.  The  visible 
white  water  far  wake  at  the  surface  usually  spreads  as  shown  at  an  angle  of  1  to  1  1/2  degrees,  down¬ 
stream  from  the  initial  spreading  region  directly  adjacent  to  the  stem  of  the  ship. 

OCEAN  WAVE  MODIFICATIONS 

Relative  to  a  body  fixed  frame,  the  far-field  surface  wake  produced  by  a  moving  vessel  is  steady 
and  can  be  approximated  as  a  parallel  shear  u  (y ) .  Here  u  is  the  velocity  component  in  the  direction 
of  ship  motion,  x ,  and  y  is  the  transverse  horizontal  coordinate.  Velocity  gradients  in  the  y -direction 
are  assumed  to  be  much  greater  than  those  in  the  x -direction.  Variations  in  the  x  direction  are 
assumed  to  be  extremely  slow.  For  surface  waves  interacting  with  a  steady  parallel  shear  current,  the 
total  or  extrinsic  frequency,  go,  and  the  component  of  wave  number  in  the  x  direction,  kx,  are  con¬ 
served  (James  et  al.,  1977).  The  random  wave  field  may  be  characterized  by  means  of  the  spectral 
density.  To  describe  the  refraction  of  these  waves  by  the  steady  surface  current,  the  frequency- 
direction  spectral  density  is  most  convenient  since  the  total  frequency  go  is  conserved.  The  frequency- 
direction  spectral  density  <I>(co,  0  ;jt)  is  defined  from 

dE  =  <f>(a ),6-,x)wdudd  (3) 

where  dE  represents  the  mean  square  surface  displacement  contribution  for  waves  propagating  in 
directions  between  0  and  6  +  dd  and  possessing  total  frequencies  between  go  and  go  +  d go.  The  total 
energy  is  obtained  by  integration  over  all  angles  —7 r  <  0  <  tt  and  over  all  frequencies  0  <  go  <  00, 
yielding 


<f2>  =  f  f  4>(co,0;x)  GodoodO. 

JwJ0 


(4) 


It  follows  that  the  contribution  to  surface  elevation  corresponding  to  frequencies  of  Aco  about  goj  , 
and  directions  A 0  about  0,  is  given  by  Wang  (1985)  as 

fy  =  V  24>(w;-  ,0,  ,x)  gOj  Ago  A0  • 


cos  -s  kjj  (x  sin  0,-  +  y  cos  0,)  —  c o;  t  +  <foj 


(5) 


where  the  wave  number  corresponding  to  refracted  waves  of  frequency  goj  traveling  in  a  direction 
6j  is  determined  from  the  dispersion  relation  for  linear  deep  water  waves  as 


kU  = 


0 )i 


V 1  +  4Ar?  -  1 
2At) 


(6) 


where 


7]  -  GO*  sin  0,-,  GO*  =  GO  Vm/g  ,  A(y)  =  u(y)/Vm 
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Fig.  5  —  A  typical  aerial  view  of  the  white  water  region  behind  a  surface  ship  underway;  from 
Peltzer  (1984). 
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and  (ftij  represents  a  random  phase  shift  uniformly  distributed  between  0  and  2it.  In  Eq.  (6)  Vm 
represents  the  maximum  value  or  strength  of  the  surface  current  distribution,  u(y),  and  g  denotes  the 
acceleration  of  gravity.  The  full  free  surface  realization  or  displacement  at  a  given  time  is 
represented  by  the  double  summation 

,  y)  =  ££  tu(a>j ,6t;x,y).  (7) 

i  j 


It  is  apparent  that  the  major  task  in  characterizing  the  interaction  of  the  ambient  wave  field  and 
the  surface  current  u(y)  is  the  evaluation  of  the  modified  spectral  density  3>.  The  local  value  of  spec¬ 
tral  density  is  determined  from  its  ratio  to  the  undisturbed  value  along  a  ray 


*  =  $00 


JL 

$ 

^  00 


=  <f>  A 

^oo  71  ? 


(8) 


where  from  conservation  of  wave  action  density 


A  = 


VTT  4Ar?  -  1 

2Ar] 


VI  +  4  At;  ' 


(9) 


The  mean  square  surface  elevation  also  can  be  determined  by  integrating  the  spectral  density  in  the 
initial  plane  (oia,,  0^)  after  a  transformation  of  variables 


<f  >  $00  A  Woo d 

Along  a  ray  Wo,  =  c o  =  constant  and 

sin  6  - 


co. 


86 

dd„ 


dd^. 


(10) 


sin  6 , 


(1  —  A  co*  sin  d,*)2 


(11) 


for  the  steady  current,  simple  shear  u(y )  case  considered  here.  For  evaluation  of  Eq.  (10),  the  spec- 

(1  -  A  r,m)5 

tral  density  ratio  A  (co^,  6X)  =  — - - - — .  The  subscript  oo  denotes  the  ambient  value  of  a 

(1  +  A  ij,*) 

quantity  along  a  ray  outside  the  region  of  current  variation. 


The  above  equations  form  the  basis  of  a  numerical  procedure  to  determine  the  surface  elevation 
field  (realization)  and  the  mean-square  surface  elevation  field.  Limits  for  the  integration  in  co  -  6 
space  are  required  to  complete  the  description.  These  limits  are  determined  from  the  conditions  of 
wave  breaking  and  wave  reflection  which  result  due  to  interaction  with  the  surface  current  pattern. 

The  wave  breaking  condition  selected  by  James  et  al.  (1977)  and  which  has  been  used  thus  far 
in  the  present  work  represents  an  extreme  condition.  The  wave  breaking  boundary  is  defined  by  the 
condition  that  the  local  component  of  the  current  in  the  direction  of  the  wave  motion  is  equal  but 
opposite  to  the  local  intrinsic  wave  group  velocity.  In  the  dimensionless  terms  of  the  present  formula¬ 
tion  this  is  equivalent  to 

7  oo  =  -™  (12) 
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or, 


V  = 


1 

4  A ' 


(13) 


These  conditions  are  conservative  and  they  are  independent  of  the  precise  value  of  the  wave  ampli¬ 
tude,  and  the  wave  amplitude  defined  by  these  conditions  is  infinite.  In  a  continuous  refraction  pro¬ 
cess,  waves  break  well  before  these  boundaries  are  reached,  and  more  energy  and  momentum  are  lost 
to  wave  breaking  than  is  indicated  by  the  above  conditions.  It  is  more  realistic  to  establish  the  wave 
breaking  boundary  in  some  manner  based  on  a  local  wave  steepness  condition  (see,  for  example, 
Ramberg  and  Griffin  (1987)  and  Ochi  and  Tsai  (1983)).  These  represent  criteria  which  depend  both 
on  the  initial  spectral  distribution  and  the  integrated  growth  of  the  particular  spectral  component. 


The  wave-breaking  boundary  is  utilized  to  establish  limits  of  integration  for  the  mean  square 
surface  elevation  energy  <p>  or,  equivalently,  limits  for  the  contributions  to  the  surface  elevation, 
Utilizing  Eq.  (13)  as  the  definition  of  the  wave  breaking  boundary  6b  the  result  is 


sin  8b  = 


1 

4Aco*  ’ 


(14) 


This  establishes  the  lower  integration  limit  for  6.  For  waves  with  components  in  the  direction  of  the 
current,  6  >  0,  interaction  with  the  current  can  produce  reflected  waves.  Waves  which  are  about  to 
reflect  are  characterized  in  the  physical  plane  at  6  =  tt /2. 


To  obtain  some  insight  into  the  spectral  distribution  throughout  the  interaction  zone  the  contribu¬ 
tions  from  various  ranges  of  wavelengths  to  the  surface  energy  must  be  quantified.  This  necessitates 
establishing  modified  limits  of  integration  for  the  energy  integral  previously  developed.  The  locus  of 
constant  wavelength  within  the  interaction  zone  in  either  the  initial  or  physical  plane  is  given  by 


sin  A,*, 


g 

A Vm  Woo 


sin  6k. 


(15) 


This  equation  provides  the  value  of  0oo  which  produces  for  a  total  frequency  a  wave  number  of  k 

at  a  position  y  where  the  local  current  is  given  by  Vm  A (y),  and  it  therefore  defines  the  required 
locus  in  the  initial  plane.  Equation  (15)  is  used  to  establish  the  required  limits  of  integration  for  the 
surface  energy  as  discussed  above.  Of  course,  the  limitations  imposed  by  Eq.  (15)  have  to  be  incor¬ 
porated  with  respect  to  those  limitations  already  associated  with  wave  breaking  and  wave  reflection. 
In  the  above  equations  represents  the  spectral  density  of  the  ambient  random  wave  field.  The 
simplest  case  corresponds  to  a  spectrum  similar  to  that  of  white-noise 


4>oo  (wo,,  0M)  =  constant, 

when  the  constant  may  be  taken  as  unity  without  any  loss  of  generality.  This  was  the  original 
approach  taken  by  James  et  al  (1977).  More  importantly,  other  random  sea  state  spectra  may  be  used 
for  such  as  those  of  Pierson-Moskowitz,  JONSWAP  etc.,  and  the  recently-proposed  model  of 
Phillips  (1985).  For  example,  the  frequency-direction  spectral  density  of  Pierson-Moskowitz  is  given 
as 


$ oo  (Wqq,  0 oo ) 


exp 


Wo 


-.74/ 


cosp 


(16) 
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where  i/w  represents  the  wind  velocity,  \pw  represents  the  wind  direction,  v  and  p  are  parameters  typ¬ 
ically  equal  to  1  and  2,  respectively,  and  C  is  a  constant.  The  co  ~6  behavior  of  Eq.  (16)  is  consistent 
with  the  spectral  density  as  defined  by  Eq.  (3).  The  various  models  in  common  use  are  described  by 
Wang  (1985)  and  Phillips  (1985). 

KELVIN  SHIP  WAVES 

In  many  surface  hydrodynamics  applications  it  is  of  interest  to  predict  the  wave-making  resis¬ 
tance  of  a  particular  ship  hull  and  the  resulting  free  surface  ship  wave  pattern.  Such  predictions  can 
be  used  for  the  analysis  and  study  of  the  hull  design  and  its  wake  characteristics.  To  meet  this  need, 
the  far-field  ship  wake  computer  program  FFSW  has  been  developed.  A  detailed  description  of  this 
program  and  discussion  of  the  theory  implemented  to  develop  the  computational  model  can  be  found 
in  Keramidas  and  Bauman  (1987). 

The  computational  model  employed  in  the  computer  program  is  based  on  Slender  Ship  theory 
(Noblesse,  1983)  with  a  number  of  improvements  over  existing  computer  programs  which  employ  the 
same  theory.  These  improvements  allow  for  a  more  accurate  description  of  hull  geometry  and  panel 
generation.  The  computer  program  is  divided  in  three  major  parts:  a)  Input,  geometry  (panel)  gen¬ 
eration  and  source  strength  generation;  b)  Wave  spectra,  wave  resistance  and  wake  calculations;  and 
c)  Output  for  printer  or  disc  file.  The  FFSW  program  is  completely  interactive  and  menu  driven,  so 
that  it  requires  only  minimum  user  interface. 

The  implementation  of  the  Slender  Ship  theory  into  a  useable  computational  tool  consists  of 
several  steps.  Each  of  these  steps  provides  the  necessary  input  to  the  computational  model  for  obtain¬ 
ing  numerical  results  for  the  wave  resistance  and  the  Kelvin  wave  field.  The  first  and  probably  the 
most  important  step  is  the  generation  of  the  surface  panels  and  the  specification  of  the  source  strength 
densities;  the  accuracy  of  the  calculations  depends  a  great  deal  on  the  accurate  representation  of  the 
ship  hull  surface  by  the  panels  as  well  as  as  the  accurate  specification  of  the  unit  normal  vector  on 
each  of  the  panels.  A  set  of  algorithms  has  been  written  to  perform  both  the  panelization  of  the  hull 
surface  and  the  calculation  of  the  unit  normal.  The  algorithm  for  the  latter  produces  results  for  both 
the  surface  area  of  each  panel  and  the  unit  normal.  This  is  a  much  improved  algorithm  over  existing 
ones,  since  it  takes  into  account  the  curvature  of  the  panel.  The  input  information  required  to  per¬ 
form  the  panelization  process  is  obtained  from  a  data  file  which  contains  the  offset  lines  of  the  partic¬ 
ular  ship  hull  under  consideration. 

Two  sets  of  calculations  are  possible;  the  first  set  of  calculations  that  one  can  perform  is  of  the 
wave  spectra,  for  a  given  ship  hull  configuration  and  for  a  specified  number  of  ship  speeds.  Wave 
making  resistance  is  then  calculated  from  the  wave  spectra  results  and  can  be  saved  as  an  output  file: 
the  second  set  of  calculations  includes  the  free  surface  elevations  of  the  far  field  wake  as  well  as  both 
x-  and  y -slopes. 

SHIP  WAVE  MODIFICATIONS 

The  theoretical  model  employed  here  is  essentially  that  of  Peregrine  (1971)  as  modified  by  Grif¬ 
fin  (1988).  It  is  assumed  that  the  relatively  short,  divergent  waves  generated  at  the  stem  of  the  ship 
(source)  move  away  from  the  center  line  of  the  wake  as  they  are  refracted  by  the  mean  flow,  and  that 
they  pass  through  only  a  small  region  of  the  wake.  The  source  is  assumed  to  represent  the  stern 
region  as  the  chief  source  of  wave  generation  in  the  pattern.  It  is  understood,  however,  that  in  many 
cases  waves  generated  at  the  ship’s  bow  also  may  contribute  to  the  overall  wave  system. 
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The  wake  is  assumed  to  be  narrow  relative  to  the  lateral  dimension  of  the  far-field  wave  pattern 
and  to  be  uniform  or  at  most  slowly-varying  with  downstream  distance  from  the  moving  source.  This 
assumption  places  stringent  limitations  on  the  allowable  extent  of  the  surface  current  profile.  In  a 
practical  sense  this  limits  the  origin  of  the  flow  to  the  wake  caused  by  the  passage  of  the  ship  (typi¬ 
cally  the  wake  is  approximately  two  to  five  ship  beams  in  width).  Naturally-occurring  surface  current 
fields  generally  are  greater  in  extent.  These  are  reasonable  assumptions  under  most  circumstances. 
The  moving  source  is  assumed  to  generate  waves  of  all  frequencies,  but  only  the  ship  waves  which 
are  steady  and  move  with  the  source  speed  are  considered  in  the  model.  It  is  the  well-known  diver¬ 
gent  cusp  waves  which  are  most  visible  near  the  boundaries  of  the  Kelvin  ship  wave  pattern,  and  it  is 
the  general  pattern  and  modifications  of  these  waves  which  are  predicted  here. 

A  coordinate  system  is  chosen  with  x  along  the  center  of  the  wake  or  the  ship  track,  y  normal 
to  the  ship  track,  and  z  vertically  downward.  The  undisturbed  free  surface  is  located  in  the  plane 
z  =  0.  The  wavemaking  source  is  assumed  to  have  a  steady  velocity  U  in  the  x  direction.  The 
wake  is  represented  in  its  most  general  form  by  a  velocity  distribution  u(y,z)  with  the  properties 
0  <u(y  ,z)<U  for  —  d  <y  <  +d  and  u  =  0  for  y  outside  of  this  region.  Assume  that  the  source¬ 
generated  waves  have  traversed  the  wake  region  and  are  propagating  away  freely  in  still  water,  i.e., 
there  is  no  appreciable  ambient  surface  wave  field.  Then  if  the  wave  pattern  is  steady  and  moving 
with  the  velocity  U ,  individual  wave  components  satisfy  the  relation 

c  =  U  cos  6 ,  (17) 

where  c  is  their  phase  velocity  and  6  is  the  angle  that  their  direction  makes  with  the  x  direction. 

The  wave  energy  is  propagated  inside  the  wake  at  a  velocity  which  is  equal  to  the  sum  of  the 
group  velocity  normal  to  the  wave  crests  and  the  local  velocity  u(y)  of  the  water.  Based  upon  the 
previous  work  of  Peregrine  (1971)  and  of  Longuet-Higgins  and  Stewart  (1961),  the  properties  of  the 
refracted  waves  inside  the  wake  are 


and 


Ci00 


c  —  u(y)  cos  6 


(18a) 


cos  0j(y)  = 


c2  COS0 


(18b) 


( c  —  u(y)  cos  6)2 

Here  the  subscript  1  denotes  values  of  the  wave  properties  inside  the  wake  and  the  absence  of  the 
subscript  denotes  the  freely  moving  waves  outside  of  the  wake. 


If  the  substitution 


W{y)  =  1  -  u(y)/U  (19) 

is  made  and  Eq.  (17)  is  introduced,  then 

Ci(y)  =  c /W(y),  cos  0](y)  =  cos  0/W2(y). 
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The  group  velocity  in  deep  water  is  equal  to  half  the  phase  velocity,  so  that  within  the  wake  a  wave 
packet  has  the  velocity 


U_  cos2  0  U_  cos  0  [lf4(y)  —  cos2  0]1/2 
2  W2(y)  ’  2  W3(y) 


If  the  source  is  at  the  origin  at  the  time  t  =  0,  then  the  waves  generated  when  the  source  was  at 
(  —  Ut ,  0  )  have  a  locus  in  the  x  ,y  plane  given  by 


x  =  -u,  +lUt  cos2  e  +  (  [(1 ^(iW  0  +  2  (1  -  wuMlrnn 

0  j  _ zi  rn/4/'  s-\  _ 2  /nl/2 


COS  0  [lT4(f)  -  COS2  0]1/2 


(20a) 


1 


y  =  d  4-  — ■  Ut  sin  6  cos  6  —  j 


W\$)  sin  6dl 


2  o  [W4tf)  -  cos2  0]1/2 

for  that  part  of  the  pattern  with  y  >  d .  The  envelope  of  the  wave  pattern  given  by  the  equations 

d 


(20b) 


=  1 


[1  -  Ik4(C)]  cos3  6  [3  cos2  6  -  2  -  2  Ik4(f)]  d? 

(2-3  cos2  0)[JP4(f)  -  cos2  0]3/2 


(21a) 


=  d  +  f  W4(f  )  sin3  0  [3  cos2  6  -  2W4]  d$ 

J  ( 2-3  cos2  0)  [W4({)  -  cos20]3/2  ’ 


(21b) 


which  are  derived  by  eliminating  the  time  t  from  Eqs.  (20)  and  taking  the  derivative  3/90  of  the 
resulting  intermediate  equation  (Peregrine,  1971).  This  envelope  is  valid  for  an  arbitrary  surface- 
plane  (z  =  0)  velocity  profile  u(y)  over  the  total  wake  width  2d.  The  refracted  steady  wave  pattern 
corresponding  to  such  a  profile  can  be  computed  numerically  using  Eqs.  (20)  to  give  the  pattern  of 
the  wave  crests.  The  corresponding  Kelvin  wave  amplitudes  for  a  single  moving  source  can  be  deter¬ 
mined  by  combining  the  approach  taken  by  Munk  et  al.  (1987)  with  the  present  one. 


Lines  of  constant  phase  in  the  wave  pattern  are  given  by  the  equation 


X  =  k  •  r  —  at  =  k(x  cos  0  +  y  sin  0  )  —  at . 
In  the  present  case  for  steady,  deep  water  waves 


c  =  U  cos  0  =  -7-  =  —  , 
k  a 

so  that  the  phase  function  is  transformed  to 


X  =  -&r  (x  cos  0  +  y  sin  0  —  ct) , 
cl 


or 


X  = 


.  6  - —  (jc  cos  0  +  y  sin  0  —  U  cos  0  t ). 
U2  cos2  0 


(22) 
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This  equation  represents  a  relation  between  t ,  6  and  X  when  x  and  y  are  substituted  from  Eqs.  (20a) 
and  (20b). 

There  are  a  number  of  features  of  the  overall  modified  ship  wave  pattern  which  can  be  obtained 
from  the  foregoing  equations.  At  the  time  of  Peregrine’s  original  work  there  were  virtually  no  exper¬ 
imental  results  with  which  to  compare  the  predictions.  This  situation  has  changed  somewhat  as  more 
wake  observations  have  become  available  from  airborne  and  spaceborne  SAR  and  other  forms  of 
imagery.  Thus  even  such  a  relatively  simple  approach  can  provide  useful  initial  information  and  pos¬ 
sible  directions  for  more  detailed  further  studies.  A  first  approximation  can  be  obtained  by  assuming 
that  the  surface  current  profile  is  given  by  u(y)  =  V9  a  constant,  for  \y  \  <d.  Then  no  numerical 
integration  is  required  and  the  equations  for  the  envelope  of  the  modified  wave  pattern  simply  reduce 
to 


=  d(  1  -  W 4)  cos3  0  (3  cos2  0  -  2  -  2 IT4) 

(2-3  cos2  0)  (W4  -  cos2  0)3/2 
and 


d  IT4  sin3  0  (3  cos2  0  -  2  If4) 

(2-3  cos2  0)  (W4  —  cos2  0)3/2  ' 


(23a) 


(23b) 


It  is  instructive  to  look  briefly  at  two  special  cases  of  these  equations.  Then  some  further  calcu¬ 
lations  can  be  made  for  some  cases  of  more  general  interest.  The  far  field  behavior  of  the  waves  is 
readily  obtained  from  the  above  equations.  The  denominator  on  the  right  hand  side  of  each  equation 
must  approach  zero  in  order  for  x  and  y  to  increase  without  bound  as  the  far  field  is  approached. 
There  are  two  ways  in  which  this  situation  can  occur. 

Case  One:  W4  >  2/3. 


The  denominators  of  Eqs.  (7a,  b)  then  are  zero  if 

(2-3  cos20)  =  0. 


(24) 


This  reduces  to  the  usual  Kelvin  ship  wake  pattern  (0  =  35.3  degrees  and  a  —  19.5  degrees)  and 
represents  what  Peregrine  has  described  as  a  ’weak’  wake.  The  condition  W4  >  2/3  implies  that 
u  <  0. If/. 


The  envelope  of  the  wave  pattern  is  bounded  by  the  asymptotes  given  by  the  line(s) 

sin  0  cos  0 


y  = 


1  +  sin2  0 


or,  in  terms  of  the  half-angle  a.  of  the  wake, 

sin  0  cos  0 


tan  a  =  — 


1  /2  tan  0 


1  +  sin2  0  1/2  +  tan2  0 


(25) 


(26) 


Here  the  angle  a  is  given  by 


a  =  arc  tan 


z 

X 
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Equation  (25)  has  been  shown  by  Crapper  (1984),  Newman  (1977)  and  by  Lighthill  (1979),  among 
others,  to  represent  the  usual  Kelvin  wave  pattern  for  a  ship  moving  in  still  water. 

Case  Two:  W4  <  2/3. 


For  this  case  the  denominator  is  zero  if 

cos2  0  =  W4 .  (27) 

This  condition  results  in  a  larger  value  of  6  than  for  the  unmodified  Kelvin  wave  pattern.  The  angle  6 
of  the  modified  divergent  wave  system  varies  over  only  a  relatively  small  range  as  compared  to  the 
previous  case  since,  at  y  =  d ,  the  limiting  value  is  given  by 

cos2  0max  =  2/3 IT4.  (28) 

This  gives  the  maximum  value  of  the  angle  6  for  a  wake  condition  represented  by  the  current  field  of 
strength  W .  The  behavior  of  the  waves  at  infinity  also  can  be  determined  in  a  similar  manner  from 
the  full  Eqs.  (20a, b)  as  shown  by  Peregrine  (1971). 

Not  only  the  angular  extent  and  phase,  but  also  the  characteristic  lengths  of  the  divergent  wave 
pattern  are  altered  by  the  presence  of  the  current.  The  modified  lengths  of  the  steady  waves  in  the 
pattern  are  given  by 


x  2t  U2  2  a 
X  =  - cos  6 . 

8 


(29) 


In  terms  of  the  reference  Kelvin  wave  length  Xq  for  a  given  set  of  conditions,  the  divergent,  small 
wavelength  portion  of  the  ship  waves  at  the  outer  edge  of  the  v  wake  is  given  by 


W4 


\)  COS20O 

where  cos2  0O  =  2/3.  Equation  (13)  can  be  recast  in  the  form 


(30) 


X/L  =  2tt  Ul  / gL 


C' 


cos20  =  27rF2cos2  6 . 


(31) 


for  a  ship  of  length  L  and  Froude  number  F  =  U  /  ^  gL  .  In  terms  of  the  wake  half-width  d. 


X  /  d  •  d  /  L  =  2-kF2  cos2  6. 

The  limiting  value  of  X  /  d  =  1,  since  within  the  limitations  of  the  present  simplified  theory  the 
wavelength  X  should  be  less  than  the  wake  half-width  d .  Then  the  above  equation  reduces  to 

d  /  L  =  2i rF2  cos2  6 . 

The  half-width  d  is  approximately  equal  to  the  ship’s  beam  B  (since,  empirically,  2d  —  2  to  3 B), 
which  gives  the  result 


B  /L  =  2irF2 cos2 d . 

For  a  typical  ship  B  /  L  =  0.10  to  0.15,  so  that  we  may  choose  as  a  representative  value 
B  /  L  =0.13.  When  the  half-angle  of  the  modified  Kelvin  wake  is  a  =  14  degrees,  the  wake  velo¬ 
city  parameter  W4  =  0.25  and  the  wave  direction  angle  is  0max  =  66  degrees.  Introducing  these 
numbers  into  Eq.  (31)  gives 
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F 2  =  0.122;  F  =  0.35 

as  the  limiting  value  of  the  Froude  number  for  these  conditions.  Thus  this  relatively  simple  theory 
can  represent  a  first  approximation  for  a  realistic  surface  ship  condition  as  represented  by  the  Froude 
number.  The  Froude  number  of  the  high-speed  destroyer  (model  and  full  scales)  is  F  =  0.28  for  the 
conditions  considered  here. 

LIMITING  WAVE  STEEPNESS 

If  the  limiting  wave  steepness  sL  is  specified  in  some  manner,  then,  for  given  ambient  condi¬ 
tions  a  km,  to oo  where  u  =  0,  the  solution  of  the  governing  model  equations  gives  the  local  surface 
shear  current  which  accompanies  the  onset  of  wave  breaking.  A  reasonable  estimate  of  the  expected 
local  steepness  s  =  ak  at  the  onset  of  breaking  is  sL  =  0.34  (Ramberg  and  Griffin,  1987).  The  clas¬ 
sical  Stokes  criterion  for  the  breaking  of  steep,  symmetric  waves  in  deep  water  has  been  computed  by 
numerous  investigators  and  is  sL  =  0.446.  This  is  another  conservative  criterion  for  predicting  the 
onset  of  breaking  since  it  represents  fewer  breaking  waves  in  a  given  record.  An  even  more  conser¬ 
vative  approach  is  that  of  James  et  al.  (1977).  The  breaking  condition  there  is  taken  to  be  the  vanish¬ 
ing  of  the  local  absolute  group  velocity  Cg .  Then  the  opposing  current  at  the  breaking  condition  is 
given  by 


u(y) sin0  =  —  cg  .  (32) 

where  cg  is  the  local  intrinsic  group  velocity  (cg  =  a / k). 

The  principal  problem  which  arises  in  trying  to  implement  any  criterion  is  the  specification  of 
the  incident  or  ambient  and  breaking  wave  steepnesses.  For  monochromatic  incident  waves  as  in  a 
laboratory  channel  the  problem  is  a  relatively  minor  one.  Both  the  initial  steepness  and  the  limiting 
steepness  at  the  onset  of  breaking  are  relatively  easy  to  measure  and  to  characterize  (Melville,  1982; 
Su  et  al.,  1982;  Ramberg  and  Griffin,  1987). 

The  onset  of  breaking  or  the  limiting  steepness  for  deep  water  waves  can  be  predicted  by  com¬ 
bining  measured  laboratory-scale  breaking  wave  limits  with  some  acceptable  form  of  the  ambient  or 
incident  spectrum.  From  a  practical  standpoint  this  approach  provides  a  reasonable  estimate  based 
upon  direct  measurements  of  individual  breaking  events.  Such  an  approach  is  based  upon  the  original 
work  of  Longuet-Higgins  (1969),  Nath  and  Ramsey  (1976)  and  of  Kjeldsen  and  Myrhaug  (1978). 
Ochi  and  Tsai  (1983)  have  measured  individual  breaking  events  for  steep  nonlinear  waves  and  have 
introduced  more  realistic  wave  spectra  into  the  formulation.  This  approach  to  specifying  the  limiting 
or  onset  condition  also  is  discussed  by  Ramberg  and  Griffin  (1987),  Huang  et  al.  (1986),  and  Srokosz 
(1986). 

Increased  steepness  and/or  breaking  of  the  waves  are  strong  contributors  to  the  enhanced  wake 
signature.  The  recent  experiments  of  Banner  and  Fooks  (1985)  have  shown  that  breaking  surface 
waves  can  significantly  influence  the  radar  cross-section  (RCS)  of  the  ocean  surface. 

A  brief  summary  of  the  approach  and  a  conceptualization  of  how  it  may  be  introduced  into  the 
present  formulation  is  given  here.  The  ratio  of  the  incident,  or  ambient,  and  breaking  wave 
steepnesses  is  given  by 


SL  aLkL 
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This  ratio  can  be  expressed  in  the  equivalent  form 


*^00  (^foo/^oo) 

IT  =  (33) 

since  the  crest-to-trough  wave  height  H  =  2a .  For  the  ambient  waves  in  deep  water  the  dispersion 
relationship  gives  =  gk  <».  At  the  onset  of  wave  breaking  the  Stokes  limiting  wavelength  is 

\L  =  (1.2g  /2-k)TI  ,  (34) 

so  that  the  ratio  of  wave  steepnesses  can  be  further  reduced  to  the  form 

-  =  - =-.  (35) 

sl  (HL/l.2T?) 

Here  Tm  and  TL  respectively  are  the  wave  periods  of  the  ambient  and  local  limiting  waves. 

The  steepness  criterion  for  the  onset  of  wave  breaking  in  deep  water  as  stated  by  Ramberg  and 
Griffin  (1987)  is 


H  >  agT2,  (36) 

where  a  is  an  empiracally  determined  threshold  parameter.  The  most  comprehensive  measurements 
available  at  the  present  time  indicate  that  a  =  0.021.  Then  the  ratio  of  wave  steepnesses  reduces  to 

Soo/sL  =  1 .2H„/agTl  (37) 

and  is  expressed  in  terms  only  of  the  threshold  parameter  a  and  a  characteristic  height  H  x  and  period 
T a,  of  the  ambient  or  incident  waves. 

Srokosz  (1986)  has  introduced  a  related  criterion  for  the  onset  of  wave  breaking,  or  for  the  lim¬ 
iting  wave  height  or  steepness.  It  is  assumed  that  the  wave  will  break  if  the  downward  acceleration 
of  the  crest  exceeds  &g ,  where  (3  is  another  positive  threshold  parameter.  Longuet-Higgins  (1963, 
1969)  has  demonstrated  on  theoretical  grounds  that  the  acceleration  at  the  crest  of  the  highest  wave  is 
0.5g,  so  that  a  reasonable  choice  based  directly  upon  theory  is  13  =  0.5.  It  is  straightforward  to  show 
that  there  is  a  correspondence  between  Eq.  (36)  and  the  maximum  acceleration  at  the  crest.  Given 
that 


H  >  agT2, 

and  assuming  again  as  a  first  approximation  that  H  =  2a ,  then  Eq.  (36)  takes  the  form 


2 a  >  ag 


C  ^  2 
I'K 

OJ 


or  equivalently  that 


a  co2  >  2i r2ag  .  (38) 

The  left-hand  side  of  this  equation  is  the  acceleration  of  the  wave  crest  to  a  good  approximation,  so 
that 


/ 3  =  27T2q:. 


(39) 
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Taking  the  experimentally-determined  value  of  a  =  0.021,  the  parameter  /3  =  0.41.  Thus  the  limit¬ 
ing  acceleration  of  the  wave  crest  at  the  onset  of  breaking  is  0.41g  based  upon  experiment,  whereas 
the  limiting  acceleration  is  0.5g  based  upon  theory  (Longuet-Higgins,  1963,  1969). 

There  are  a  number  of  ways  to  specify  H „  and  Tx.  Ochi  and  Tsai  (1983)  set  Hx  proportional 
to  the  rms  wave  height  and  Tx  equal  to  the  average  time  between  positive  maxima  for  a  narrow-band 
spectrum.  Kjeldsen  and  Myrhaug  (1978)  simply  set  Hx  and  Tx  equal  to  the  rms  wave  height  and 
period  of  the  wave  field.  Unlike  the  simple  criterion  of  Eqs.  (12)  and  (13),  the  effects  of  wave 
breaking  are  not  confined  to  a  local  region  in  wave  number  or  frequency-direction  space.  At  present 
we  are  investigating  means  to  introduce  this  type  of  non-local  criterion  into  the  wave-current  interac¬ 
tion  model,  and  a  brief  discussion  is  given  later  in  the  report. 

DISCUSSION  OF  RESULTS 

Kelvin  Ship  Waves 

The  results  to  be  presented  in  this  section  include  a  brief  discussion  of  the  far  field  Kelvin  waves 
generated  by  a  ship  in  forward  motion.  The  formulation  used  for  calculating  the  Kelvin  wave  field  is 
based  on  Slender  Ship  theory  (Noblesse,  1983)  and  numerical  results  were  obtained  using  the  com¬ 
puter  program  FFSW  (Keramidas  and  Bauman,  1987).  The  numerical  results  from  FFSW  provide 
predictions  for  three  quantities  relevant  to  the  far  field  Kelvin  wave  field;  these  are  the  wave  eleva¬ 
tions,  and  the  x-  and  y -slopes  of  the  wave  elevations. 

The  availability  of  recent  model-scale  experimental  data  for  the  high  speed  destroyer  hull  with  a 
sonar  dome  (Model  5415)  makes  it  opportune  to  use  this  generic  hull  in  the  present  study  to  predict 
its  Kelvin  wave  field  and  its  wave  modifications.  The  first  set  of  results  that  is  presented  in  this 
report  was  obtained  from  the  free  surface  calculation  using  FFSW.  In  Fig.  6  the  results  obtained  for 
the  far  field  Kelvin  waves  are  shown.  The  calculations  extend  downstream  from  the  stern  of  the  ship 
to  about  five  ship  lengths  and  two  ship  lengths  on  each  side  in  the  cross  stream  direction.  The  results 
were  obtained  for  a  full-scale  speed  of  20  knots  (10  m/s),  which  corresponds  to  the  model  experiment 
(see  Table  1).  Figure  6a  shows  the  actual  wave  elevations  in  a  perspective  view,  and  Fig.  6b  shows 
the  contours  of  constant  wave  amplitude.  The  second  set  of  results  is  given  in  terms  of  wave  slopes. 
As  was  mentioned  above,  the  program  output  generates  the  x  and  y  slopes  of  the  wave  elevations, 
which  are  proportional  to  the  respective  directional  wave  numbers.  In  Fig.  7  the  total  slopes  are 
presented,  or  the  magnitudes  of  the  product  of  the  wave  number  and  the  amplitude.  These  results  are 
more  relevant  to  radar  imagery  than  are  the  actual  wave  elevations.  By  comparing  Figs.  6b  and  7 
one  can  see  that  they  have  different  characteristics  and  that  certain  features  of  the  Kelvin  wave  field 
are  present  only  in  Fig.  7. 

From  recent  Kelvin  wave  observations  the  presence  of  the  narrow  v  angle,  smaller  than  the 
usual  19.5  degree  value,  has  been  observed  as  described  earlier  in  the  report.  In  Fig.  7  one  can 
notice  the  presence  of  such  an  angle  being  formed  by  the  contours  of  the  wave  slopes  between  ak  = 
0.14  and  0.34.  This  result  is  encouraging  and  may  have  some  bearing  on  the  overall  v  wake  ques¬ 
tion,  but  further  investigations  are  required.  Another  possible  cause  of  the  narrow-v  angles  observed 
in  the  radar  imagery  is  discussed  in  the  next  section  of  the  report.  The  final  set  of  results  for  the  Kel¬ 
vin  wave  field  of  the  high  speed  destroyer  is  presented  in  Fig.  8  in  terms  of  two  wave  cuts;  the  first 
about  one  beam  distant  from  the  centerline  of  the  wake  and  the  second  at  two  beams  distant.  These 
results  illustrate  the  magnitude  of  the  wave  field  generated  by  the  passage  of  the  ship,  at  the  forward 
speed  of  20  kt  (10  m/s). 
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X/L  =  1  (STERN),  y/L  =  2 


Fig.  6(a)  —  The  Kelvin  ship  wave  pattern  in  the  wake  of  a  twin-screw  high  speed  destroyer  in 
steady  motion  at  a  forward  speed  of  20  kt  (10  m/s). 

A  three  dimensional  surface  plot  of  the  wave  field. 
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Fig.  6(b)  —  The  Kelvin  ship  wave  pattern  in  the  wake  of  the  a  twin-screw  high  speed  destroyer  in 
steady  motion  at  a  forward  speed  of  20  Kt  (10  m/s). 

Contours  of  constant  wave  amplitude.  (1  ft  =  0.305  m) 
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Fig.  7  —  A  contour  plot  of  the  wave  slopes  ak  in  the  wake  of  a  high  speed  destroyer  at  20  kt  (10 
m/s).  (1  ft  =  0.305  m) 


Wave  Rmpli'tude,  C-F*t) 


Fig.  8  —  Longitudinal  wave  cuts  in  the  Kelvin  wake  of  a  high  speed  destroyer  at  20  kt  (10  m/s).  (1 
ft  =  0.305  m) 
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Table  2  —  Ship  Wave  Modification  by  a  Surface  Wake  Flow  Field 

(Divergent  Waves) 


Wake  Half-angle,  a 
(deg) 

Relative  Current,  V/U 

Ship  Speed,  V 
(kt) 

Wavelength,  X 
(m) 

16 

0.23  (W4  =  0.35)+ 

12 

8.48  (X/Xo  =  0.53) 

8 

3.76 

4 

0.95 

14 

0.29  (W4  =  0.25) 

12 

6.20  (X/X o  =  0.38) 

8 

2.75 

4 

0.69 

12 

0.37  (W4  =  0.16) 

12 

3.84  (X/X0  =  0.24) 

8 

1.71 

4 

0.43 

6 

0.55  (IF4  =  0.043) 

12 

1.02  (X/X0  =  0.064) 

8 

0.46 

4 

0.11 

4 

0.63  (IF4  =  0.029) 

12 

0.47  (X/X0  =  0.029) 

8 

0.21 

4 

0.052 

+  W  =  1  -  V/U 


I 

I 

I 
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Fig.  10  —  Typical  wavelengths  of  the  divergent  waves  in  the  modified  surface  ship  wake  as  a  func¬ 
tion  of  the  current/ship  speed  ratio  V/U.  For  V/U  <  1/10,  the  wavelengths  are  those  of  the  unmo¬ 
dified  Kelvin  wake. 
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NUMBER  OF  OBSERVATIONS 


15 


- SEAS  AT  ( ~  1 978) 

- SIR-B  (1984) 
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Fig.  11  —  Observed  distribution  of  SEASAT  and  SIR-B  ship  wake  angles  (SEAS AT  results  adapted 
from  Hammond  et  al.,  1985).  The  SIR-B  results  were  provided  by  Dr.  A.  M.  Reed  of  DTRCEN 
(1986). 
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Ambient  Wave  Modifications 


The  results  to  be  discussed  in  this  section  of  the  report  include  sample  calculations  of  the  modif¬ 
ication  of  a  random  deep  water  wave  field  by  the  presence  of  a  surface  current  pattern  U(y).  The 
general  features  and  dimensions  of  the  surface  region  containing  the  current  correspond  to  the  surface 
wake  of  a  full-scale  destroyer.  The  width  of  the  current  field  is  56  meters.  This  is  nominally  2.5 
ship  beams,  for  B/L  =  0.133  as  in  the  case  of  the  high  speed  destroyer  of  length  L  =  143  m  (470 
ft).  The  maximum  value  of  the  surface  current  is  chosen  to  be  VM  =  75  cm/sec  (1.5  kt),  which  is 
representative  of  the  full-scale  ship  wake  as  scaled  up  from  the  model  tests  and  computations.  This  is 
7.5  percent  of  the  full-scale  ship  speed  of  20  kt,  a  condition  which  corresponds  to  the  surface  mean 
flow  (relative  to  still  water)  at  one-half  ship  length  astern.  The  results  largely  have  been  chosen  to 
demonstrate  the  effectiveness  of  the  wave-current  interaction  in  eliminating  the  short  wave  com¬ 
ponents.  For  example,  the  wave  number  range  of  k  =  0.2  to  1  cm-1  in  deep  water  is  comparable  to 
the  operating  wave  number  ranges  of  L-,  C-  and  A -band  radars. 

A  Sea  State  Four/Five  ambient  wave  field  was  chosen  for  the  computations  using  a  Pierson- 
Moskowitz  spectrum  with  cosine-squared  spreading  as  an  example.  The  computational  grid  in  all  of 
the  examples  consisted  of  an  array  of  51  X  51  points.  The  nominal  wind  speed  used  in  the  simula¬ 
tion  was  U  —  12  m/sec  (24  kt).  The  patch  of  ocean  surface  measured  70  m  x  70  m,  with  the  sur¬ 
face  wake  or  current  field  covering  the  central  56  m  in  the  y  or  cross-wake  direction.  This 
corresponds  nominally  to  the  wake  dimensions  discussed  in  the  previous  paragraph,  although  the 
simulation  is  somewhat  of  an  idealization  because  a  simple  parabolic  variation  of  the  surface  mean 
velocity  over  the  surface  wake  was  employed.  The  horizontal  x  and  y  directions  are  scaled  by  the 
factor  g/Vh- 

Figure  12  shows  the  surface  wave  pattern  for  wave  numbers  between  k  =  0.2  and  1cm-1  when 
the  prevailing  wave  system  is  in  the  same  direction  as  the  surface  flow  pattern.  This  is  the  most 
effective  relative  orientation  of  the  wave  and  current  systems  for  reducing  the  short  wavelengths 
(Cooper  et  al.,  1987).  The  modified  wave  pattern  includes  all  wave  numbers  greater  than  k  =  0. 
The  maximum  peak-to-trough  wave  height  /7rms  =  1.1  cm  is  over  the  limited  short  wave  number 
range  investigated  in  this  example.  There  is  a  clear  reduction  in  the  short  wavelengths  over  the  extent 
of  the  surface  wake  pattern  as  compared  to  the  region  outside  of  the  wake. 

When  the  wave  number  interval  is  expanded  to  k  =  0.0005  to  3  cm-1,  the  surface  wave  pattern 
is  as  shown  in  Fig.  13.  The  surface  corresponding  to  all  wave  numbers  greater  than  k  =  0  now 
takes  on  the  appearance  of  a  random  sea,  with  a  maximum  wave  height  Hm s  =  121.2  cm.  There  are 
no  visually  identifiable  regions  of  wave  modification  in  the  surface  pattern,  as  in  the  short  wavelength 
case,  but  the  wave  spectra  in  the  current  field  are  modified  as  shown  below. 

The  strength  of  the  surface  wake  diminishes  sharply  as  the  distance  astern  of  the  ship  is 
increased.  At  a  downstream  distance  of  approximately  twelve  ship  beams  (1.5  ship  lengths),  the  max¬ 
imum  surface  mean  velocity  is  about  one  percent  of  the  ship  speed,  or  10  cm/sec  (0.2  kt).  The  sur¬ 
face  pattern  in  the  short  wave  regime  (k  =  0.2  to  1  cm-1)  is  shown  in  Fig.  14,  again  for  all 
wavelengths  in  the  modified  pattern  greater  than  k  =  0.  There  is  little  evidence  of  a  diminished  short 
wave  population  in  the  wake  region,  so  that  a  reduction  in  the  radar  backscatter  from  the  surface  is 
not  likely  to  occur  under  these  conditions. 

Another  example  is  shown  in  Fig.  15  for  the  case  when  the  directions  of  the  surface  current 
field  and  the  wave  system  are  oriented  at  ninety  degrees  to  each  other.  In  this  case  there  is  no 
appearance  of  modification  in  the  short  waves  in  the  range  k  =  0.2  to  1  cm-1,  although  there  is 
more  energy  in  the  system  as  denoted  by  the  maximum  waveheight  in  the  field  =  2.94  cm. 
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WAVE  DIRECTION 


Hrms  =  1.09  cm  (MAX) 


Fig.  12  —  Spatial  realization  of  the  modified  wave  elevation  in  a  surface  wake  or  current  field.  The 
maximum  surface  velocity  is  VM  =  75  cm/sec.  The  horizontal  field  is  70  m  in  each  direction,  with 
the  current  field  extending  over  the  central  56  m.  The  ambient  spectral  density  is  Pierson- 
Moskowitz,  with  a  wind  speed  of  12  m/sec  (24  kt).  All  wave  numbers  greater  than  k  =  0  in  the 
modified  field  are  shown,  for  incident  wave  numbers  between  k  =  0.2  and  1  cm-1.  The  waves  and 
current  are  collinear  and  unidirectional. 


PREDOMINANT 
WAVE  DIRECTION 


Hrms  =  121.2  cm  (MAX) 


CURRENT 

DIRECTION 


70  m  (x  AND  y  DIRECTIONS) 


Fig.  13  —  Spatial  realization  of  the  modified  wave  elevation  in  a  surface  wake  or  current  field. 
Except  for  incident  wave  numbers  between  k  =  0.0005  and  3  cm-1,  the  conditions  are  the  same  as  in 
Fig.  12. 
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Fig.  14  —  Spatial  realization  of  the  modified  wave  elevations  in  a  surface  wake  or  current  field. 
Except  for  a  maximum  surface  velocity  of  VM  =  10  cm/sec,  the  conditions  are  the  same  as  in  Fig. 
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Fig.  15  —  Spatial  realization  of  the  modified  wave  elevations  in  a  surface  wake  or  current  field.  The 
principal  direction  of  the  wave  field  is  normal  to  the  current;  other  conditions  are  the  same  as  in  Fig. 


12. 
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This  may  be  compared  to  the  result  in  Fig.  12  when  the  waves  and  current  were  collinear  and  uni¬ 
directional. 

Some  sample  calculations  now  are  given  for  the  frequency-direction  spectral  density  in  the  modi¬ 
fied  surface  wave  pattern.  Note  that  the  form  of  the  spectral  density  follows  that  of  Eqs.  (3)  and 
(16),  i.e.  the  frequency  varies  as  co-6.  The  direction  varies  from  —  7t/2  to  +  7r/2,  and  the  radian 
frequency  co  is  scaled  by  VM/g .  The  frequency  and  direction  are  discretized  into  25  and  24  points, 
respectively.  Figure  16  qualitatively  shows  the  variation  in  the  spectral  density  across  the  wake  when 
the  wave  system  and  the  surface  wake  are  collinear  and  unidirectional.  The  range  of  wave  numbers 
in  the  initial  plane  is  k  =  0.0005  to  3  cm-1,  and  the  modified  spectra  correspond  to  all  wave 
numbers  greater  than  k  =  0.  The  peak  value  in  the  spectral  density  is  smallest  at  the  wake  center 
line  y  =  0  (the  location  of  the  maximum  current)  and  increases  as  the  cross-wake  displacement  is 
increased.  The  surface  plot  in  Fig.  13  gives  little  evidence  of  the  wave-current  interaction,  but  the 
spectral  densities  indicate  the  modification  of  the  wave  field  as  the  wake  region  is  traversed  in  this, 
the  strongest  wake- wave  interaction.  The  distribution  of  the  wave  energy  is  confined  to  angles 
between  0  and  +  7 r  with  the  predominant  direction  of  the  incident  waves  along  the  direction 
6  =  +  7t/2. 

When  the  directions  of  the  wave  and  current  fields  are  normal  to  one  another,  there  is  little 
interaction  between  the  two.  This  is  demonstrated  by  the  two  spectra  in  Fig.  17,  one  from  the  center 
line  of  the  wake  and  one  displaced  in  the  cross-wake  direction.  There  is  only  a  slight  difference  in 
the  maxima  in  the  spectra,  but  the  angular  distribution  of  the  modified  wave  pattern  is  greater  than  in 
the  previous  example.  The  wave  pattern  is  primarily  in  the  direction  8  =  0  when  the  waves  are  nor¬ 
mal  to  the  wake  flow. 

Steep  and  Breaking  Waves 

The  relatively  simple  wave  breaking  criterion  represented  by  Eqs.  (12),  (13)  and  (32)  yields  a 
conservative  estimate  of  the  number  of  waves  in  a  given  record  which  are  at  the  onset  of  breaking, 
since  physically  this  onset  condition  occurs  before  the  intrinsic  group  velocity  c  becomes  equal  and 
opposite  to  the  current  component  in  the  wave  direction  u(y )  sin  8  However,  a  dilemma  arises  in 
attempting  to  introduce  a  steepness  criterion  such  as  Eq.  (36)  into  a  local  region  of  wave  number  or 
frequency-direction  space.  This  is  because  the  wave  height  H,  and/or  the  wave  period  T  are  averaged 
globally  in  some  manner  as  discussed  earlier  in  the  report. 

One  possible  approach  is  to  define  the  wave  steepness  in  a  local  spectral  representation  (R.  A. 
Skop,  1987).  This  can  be  accomplished  in  the  following  way.  If  the  directional  wave  number  spectral 
energy  density  is  given  by  S  (k,  6),  then  the  local  wave  steepness  is  given  to  a  good  approximation  by 

kP  S(k,  8)  ~  a2k2.  (40) 

where  p  is  an  as-yet  undetermined  integer.  When  the  onset  condition  for  wave  breaking,  H  =  agT2 
or  ak  =  0.34,  is  introduced,  the  latter  expression  reduces  to 

kp  S(k ,  6)  ~  0.116, 

which  can  be  applied  locally  in  wave  number-direction  space.  A  related  formulation  can  be  intro¬ 
duced  in  frequency-direction  space.  The  wave  number  k  is  appropriately  interpreted  as  the  resultant 
value  in  the  local  wave  direction.  Steps  are  being  taken  at  the  present  time  to  introduce  such  a  less 
conservative  onset  criterion  for  wave  breaking  into  the  wave- wake  interaction  model. 
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Fig.  16a  —  Frequency-direction  spectral  density  in  the  modified  wave  field.  The  maximum  surface 
velocity  is  VM  =  75  cm/sec.  The  ambient  spectral  density  is  Pierson-Moskowitz,  with  a  wind  speed 
of  12  m/sec  (24  kt).  All  wave  numbers  greater  than  k  =  0  in  the  modified  spectrum  are  shown,  for 
incident  wave  numbers  between  k  =  0.0005  and  3  cm-1.  The  waves  and  current  are  collinear  and 
unidirectional. 

y  =  -24  m 
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Fig.  16b  — ■  Frequency-direction  spectral  density  in  the  modified  wave  field.  The  maximum  surface 
velocity  is  =  75  cm/sec.  The  ambient  spectral  density  is  Pierson-Moskowitz,  with  a  wind  speed 
of  12  m/sec  (24  kt).  All  wave  numbers  greater  than  k  =  0  in  the  modified  spectrum  are  shown,  for 
incident  wave  numbers  between  k  =  0.0005  and  3  cm-1.  The  waves  and  current  are  collinear  and 
unidirectional. 
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Fig.  16c  frequency-direction  spectral  density  in  the  modified  wave  field.  The  maximum  surface 
velocity  is  VM  =  75  cm/sec.  The  ambient  spectral  density  is  Pierson-Moskowitz,  with  a  wind  speed 
of  12  m/sec  (24  kt).  All  wave  numbers  greater  than  k  =  0  in  the  modified  spectrum  are  shown,  for 
incident  wave  numbers  between  k  —  0.0005  and  3  cm  The  waves  and  current  are  collinear  and 
unidirectional. 
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Fig.  17a  —  Frequency 7direction  spectral  density  m  the  modified  wave  field.  The  principal  direction 

the  wave  field  is  normal  to  the  current;  other  conditions  are  the  same  as  in  Fig.  16. 
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Fig.  17b  —  Frequency-direction  spectral  density  in  the  modified  wave  field.  The  principal  direction  of 
the  wave  field  is  normal  to  the  current;  other  conditions  are  the  same  as  in  Fig.  16. 
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SUMMARY 


Remote  synthetic  aperture  radar  (SAR)  images  of  surface  ship  wakes  often  show  that  the  usual 
wedge-shaped  steady  Kelvin  ship  wave  pattern  with  a  half-angle  of  19.5  degrees  is  modified  such  that 
a  bright  narrow-v  wake  image  is  formed  downstream  of  the  ship.  The  half-angle  of  this  narrow  wake 
has  been  observed  in  the  SEAS  AT  and  DREP  experiments  to  be  as  small  as  3  to  4  degrees,  and  half¬ 
angles  of  14  to  16  degrees  are  common  to  the  SIR-B  experiments.  These  modified  Kelvin  wakes  are 
usually  observed  when  the  ship  is  traveling  in  the  azimuthal  direction  of  up  to  45  degrees  relative  to 
the  aircraft,  spacecraft  or  satellite  which  is  carrying  the  radar;  that  is,  in  a  direction  nearly  parallel  to 
the  flight  direction.  Then  the  radar  is  sensitive  to  waves  which  travel  in  the  range  direction,  i.e. 
toward  or  away  from  the  ship’s  track.  These  are  thought  to  be  ship-generated  Bragg  waves  for  the 
very  narrow  (3  to  4  degree)  wakes. 

In  contrast  to  the  bright-v  wake  images  there  also  is  a  relatively  narrow,  stem-like  dark  image 
along  the  ship  track  which  often  persists  for  many  kilometers  downstream.  The  dark  portion  of  the 
radar  image  suggests  that  this  region  is  relatively  smooth  and  free  of  short- wavelength  scatterers. 
This  region  often  is  described  as  the  turbulent  wake  or  scar  of  the  ship,  and  the  lack  of  scatterers  is 
ascribed  to  the  interaction  of  the  ship’s  turbulent  wake  with  the  surface  wave  field.  Another  sugges¬ 
tion  is  that  the  persistence  of  the  dark  scar  region  is  due  to  the  redistribution  of  surfactant  materials  at 
and  near  the  water  surface  as  a  result  of  the  ship’s  passage. 

A  model  has  been  outlined  in  this  paper  for  characterizing  the  modification  of  a  random  ocean 
surface  wave  field  by  a  surface  current  or  wake  flow  pattern.  At  the  most  fundamental  level,  for  a 
white-noise  ambient  spectrum,  the  mean-square  surface  elevation  energy  can  be  computed  for  a  wide 
range  of  conditions.  More  realistic  ambient  spectral  distributions  such  as  Pierson-Moskowitz  and  that 
of  Phillips  (1985)  also  can  be  introduced,  the  results  can  be  subdivided  into  various  wavelength  and 
frequency  ranges,  and  changes  in  the  surface  elevations  and  wave  spectra  throughout  the  wave- wake 
interaction  region  can  be  determined. 

An  important  component  in  the  overall  wave-current  interaction  model  is  the  definition  of  a  con¬ 
dition  or  criterion  for  the  limiting  steepness  or  the  onset  of  wave  breaking.  A  new  method  for 
estimating  the  limiting  steepness  from  the  recent  experiments  of  Ramberg  and  Griffin  (1987)  and  of 
Ochi  and  Tsai  (1983)  is  described  here.  This  limiting  steepness  condition  has  been  employed  recently 
by  Huang  et  al.  (1986)  in  developing  an  analytical  model  for  oceanic  whitecap  coverage. 

A  model  for  predicting  the  modification  of  the  Kelvin  ship  wave  pattern  by  interaction  with  a 
surface  current  field  also  is  described  briefly.  Not  only  are  the  lines  of  constant  phase  and  the  bounds 
of  the  Kelvin-v  angle  modified,  but  also  the  wave  lengths  and,  potentially,  the  wave  steepness  in  the 
ship  wave  pattern.  The  surface  wake  or  current  conditions  required  to  produce  a  modified  v  angle  of 
16  degrees  are  qualitatively  in  agreement  with  recent  measurements  of  the  surface  mean  flow  pattern 
in  the  near  wake  of  a  high  speed  destroyer  model. 

Recent  computations  by  Swean  (1987)  to  simulate  the  turbulent  wake  of  a  high  speed  destroyer 
from  model  test  data  have  been  combined  here  with  the  surface  wave  modification  model  to  predict 
the  interaction  of  short  waves  with  the  surface  mean  wake  flow  pattern.  The  short  waves  in  the  5  to 
30  cm  wavelength  band  interact  strongly  with  the  wake  close  to  the  ship  stern  where  the  surface  mean 
flow  is  approximately  75  cm/sec  (1.5  kt)  full-scale,  but  the  interaction  is  greatly  diminished  when  the 
surface  flow  decreases  to  about  10  cm/sec.  The  latter  condition  corresponds  to  a  distance  of  only 
about  12  ship  beams,  or  1.5  ship  lengths,  astern.  Thus  the  downstream  extent  of  the  mean  turbulent 
surface  wake  caused  by  the  ship  passage  is  not  sufficient  to  explain  the  persistence  of  the  dark  surface 
scar  for  the  observed  extent  of  often  many  kilometers. 
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APPENDIX 


The  FSWCI  Computer  Code 

The  FSWCI  computer  code  has  been  developed  to  provide  a  computation  of  the  modified  surface 
amplitude  field  in  deep  water  that  results  from  the  interaction  between  a  surface  shear  current  or  wake 
flow  pattern  and  an  ambient  wave  field.  The  basic  formulation  was  given  by  James  et  al.  (1977)  and 
is  based  upon  conservation  of  the  wave  action,  or  the  ratio  of  the  wave  energy  density  and  the  local, 
or  intrinsic,  wave  frequency.  The  original  formulation  of  James  et  al.  was  based  upon  a  white  noise 
or  “top  hat”  constant  ambient  wave  spectral  density;  in  the  present  FSWCI  code  the  white  noise 
spectrum  is  replaced  by  a  more  general  Pierson-Moskowitz  spectral  density.  Other  spectra  in  com¬ 
mon  usage  can  be  introduced  into  the  model  with  equal  facility.  The  present  version  of  the  program 
is  based  on  a  frequency -direction  formulation.  A  wave  number-direction  formulation  including  radar 
cross-section  (RCS)  predictions  in  the  perturbed  surface  plane  is  presently  under  development  for 
implementation  into  a  related  computer  code  (Skop  and  Leipold,  1988). 

The  modified  surface  wave  field  can  be  computed  over  an  arbitrary  range  of  ambient  wave 
numbers  using  up  to  51  x  51  grid  points  in  the  surface  (x  ,y)  plane.  The  local  modified  frequency- 
direction  spectral  density  can  be  computed  using  up  to  25  X  24  grid  points,  respectively.  The  pro¬ 
gram  is  menu-driven  and  user-friendly,  and  is  configured  to  run  on  the  Hewlett-Packard  9000  Series 
of  computers  (presently  Models  550  and  320).  A  listing  of  the  FSWCI  computer  code  follows. 
Further  information  is  available  from  the  first  author. 
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**  ** 

*****************************************. ******************************* 

**  ** 

PROGRAM  FSWCI 

**  ** 

************************************************************************ 
**  ** 

COMMON  /INP/  GEE , UM, YO , Y1 , W1 , PI , YM, WD, DY , WL, XS ( 100 ) ,YS(100) , 

U(100) , NX , NY , PSIW, UW , NW, NT, DELTH , AKL, NREC , NS ( 10 ) , 
TITLE, SPEC (10) 

COMMON  /OUT1/  DELW (99 ) ,THET(99) , THETA (99 ) ,WS(99) ,THIR(99) ,THL(99) 
COMMON  /OUT2/  ZETA ( 101 , 101 ) , ZETA2 ( 101 , 101) , XX ( 101 ) , YY ( 101 ) 

COMMON  /OUT3/  F (99 , 99 ) , ETA (99 , 99 ) , AK (99 , 99 ) , PH (99 , 99 ) , F2 (99 , 99 ) 


** 


** 


CHARACTER* 1  ANS 
CHARACTER* 4  0  SPEC 
CHARACTER* 80  TITLE 


10  CALL  INPUT 

*  * 

***** 

** 

CALL  SPECMOD 


Compute  the  spectral  density  at  each  point  along  the  Y-axis. 


** 


** 


** 


CALL  CLEAR_DI  SPLAY 
WRITE ( * , ' ( 5 (/ ) ,5X, 

•  u ********************************************************  ii 

./,5X,"**  Would  you  like  to  run  another  case?  (y,n):  ",$)') 

READ (* , 1 (Al) 1 )  ANS 

IF  (ANS  .EQ.  'Y'  .OR.  ANS  . EQ .  'y')  GOTO  10 

STOP 

END 


************************************************************************ 


SUBROUTINE  INPUT 


** 


** 


** 
*  * 
*  * 
*  * 


************************************************************************ 
** 

**  Dimensional  form  of  input,  defines  KL,  WL,  etc., 

**  to  be  used  with  REFRACT  and  SPECMOD. 

** 

**  INPUT  requires  dimensional  values  of  Kl,  K2 ,  and  KL,  along  with  ** 
**  dimensional  value  of  UM  in  cm/sec,  UW  in  knots,  and  all  K  ** 

**  values  in  cm-1.  ** 

**  ** 

************************************************************************ 
**  ** 

COMMON  /INP/  GEE , UM, Y0 , Y1 , W1 , PI , YM , WD, DY , WL, XS ( 100) ,YS(100) , 
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U(100) , NX , NY , PSIW , UW , NW , NT , DELTH , AKL, NREC , NS ( 10 ) , 
TITLE, SPEC(IO) 

* 

CHARACTER* 1  ANS 

CHARACTER* 4 0  NFILE , SPEC , NSPEC 

CHARACTER* 8 0  TITLE, NTITLE 

* 

PI  =  2 . 0*ASIN (1.0) 

GEE  =  980.0 

* 

CALL  CLEAR_DISPLAY 
WRITE (*,  1  ( 5  (/ )  ,5X, 

. / , 5X, " **  PROGRAM  REFRACT  **", 

. / , 5X, "********************************************************") 1 ) 
CALL  PAUSE_SEC ( 5 ) 

* 

CALL  CLEAR_DISPLAY 
WRITE (*, 1 (////, 

./, 5X, "********************************************************" f 
./f5X,"**  Would  you  like  to  read  an  input  file?  (y/n) :  ",$)') 

READ ( * , ' (Al) ' )  ANS 

IF  (ANS  .NE.  'Y'  .AND.  ANS  . NE .  'y')  GOTO  10 

* 

5  WRITE ( * , ' ( 5X, 

./f5X," **  Enter  the  name  of  the  input  file:  ",$)') 

READ ( * , ' (A40) ' )  NFILE 

OPEN (UNIT=40,FILE=NFILE,FORM=' UNFORMATTED' , STATUS= ' OLD ' ,ERR=5) 

READ (40)  TITLE , AKL, AK1 , AK2 , UW, PSIW , UM , U0 , Y0 , Y1 , NX , NY , NY0 , NY1 , 

NW, NT, NREC, (NS (N) ,N=1,NREC) , (SPEC(NN) ,NN=1,NREC) 

CLOSE (UNIT  =40) 

WRITE (*, ' (5X, 

*  f 

. / , 5X , "**  Would  you  like  to  view  or  change  any  of  the  data?  ",$)') 
READ ( * , ' (Al) ' )  ANS 

IF  (ANS  .NE.  ' Y '  .AND.  ANS  . NE  .'y')  GOTO  20 

* 

10  CALL  CLEAR_ DISPLAY 
WRITE (*, 1010)  TITLE 
1010  FORMAT (////, 

. / , 5X, " *********************************************************  t 
. /, 5X, "**  Enter  the  title  of  the  output  files.  TITLE  =  **", 

• / , A80 , 

./ , 5X, "**  ENTER  TITLE:  ",$) 

READ ( * ,  1  ( A8  0 )  ' )  NTITLE 

IF  (NTITLE  .NE.  '/')  TITLE  =  NTITLE 

* 

WRITE(*, 1020)  AKL 
1020  FORMAT (////, 
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. / , 5X,  '•  ******************************************************** 
./,5X,"**  Enter  the  minimum  wave  number  in  cm-1  for  files  **", 
. / , 5X ,  "**  Specden,  Surf 2 ,  and  Surface.  AKL  =  ",F10.5,3X, 

. /, 5X, "**  ENTER  AKL:  "  ,$) 

READ ( * , * )  AKL 


** 


WRITE (*, 1030)  AK1 ,  AK2 
1030  FORMAT (////, 

./,5X, "********************************************************", 


/,5X, 

»»  ** 

Enter 

the 

minimum  and  the 

maximum  wave  numbers 

in  **M, 

/,5X, 

n  *  * 

cm-1 

for 

the  free  stream 

spectrum  (ambient) . 

*  *  «> 

i 

/,5X, 

»  ** 

AK1 , 

AK2 

=  " , 2F10. 5, 18X, 

»  **  >» 

r 

/  /  5X , 

»  ** 

ENTER 

AK1 

,  AK2:  ",$) 

READ ( * , * )  AK1 , AK2 

U KNOTS  =  UW/51.4 
WRITE (*, 1040)  U KNOTS 
1040  FORMAT (////, 

./,5X," ******************************************* *************" 
./,5X,"**  Enter  wind  velocity  in  knots.  UW  =  ",F10.5,4X,  "**" 

./, 5X, "**  ENTER  UW:  ",$) 

READ ( * , * )  U KNOTS 
UW  =  UKN0TS*51 . 4 


PANGLE  =  PSIW* 180. 0/PI 
WRITE (*,1050)  PANGLE 
1050  FORMAT (////, 

./ , 5X, " ******************************************************** "? 
./,5X,"**  Enter  the  wind  direction  in  degrees.  **", 

. / , 5X, " **  PSIW  =  " , F10. 5, 32X,  »**", 

./ , 5X, "**  ENTER  PSIW:  ",$) 

READ ( * , * )  PANGLE 
PSIW  =  PANGLE*PI/180 . 0 

c* 

WRITE(*, 1060)  UM 
1060  FORMAT (////, 

./,5X," ********************************************************", 
. / , 5X, "**  Enter  the  maximum  current  velocity  in  cm/sec  **", 

. / , 5X, " **  (parabolic).  UM  =  ",F10.5,20Xf  "**", 

./ , 5X, "**  ENTER  UM:  ",$) 

READ ( * , * )  UM 

t  * 

WRITE (*, 1070)  U0 
1070  FORMAT (////, 

. / , 5X, " ******************************************************** »; 
./ , 5X, "**  Enter  the  velocity  of  the  constant  stream  in  **", 

./ ,  5X,  "**  cm/sec.  U0  =  ",F10.5,25X, 

./ , 5X, "**  ENTER  U0:  ",$) 

READ ( * , * )  U0 
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WRITE ( * , 1080 )  Y0,  Y1 
1080  FORMAT (////, 

./,5X, "******************************************************** u 
./,5X,"**  Enter  the  boundaries  of  the  complete  stream  in  **" 
. / ,  5X ,  "**  meters.  Y0,  Y1  =  "  ,  2F10 . 5 , 11X,  "**» 

./, 5X, "**  ENTER  Y0,  Yl:  ",$) 

READ ( * , * )  Y0 , Yl 


** 


WRITE (*,1090)  NX,  NY 
1090  FORMAT (////, 

. / , 5X, " ********************************************************* 

• / , 5X , " **  Enter  the  number  of  points  in  the  X  and  the  Y  **"^ 

./ , 5X, "**  directions.  NX,  NY  =  ",2I5,17X,  "**", 

./ , 5X, "**  ENTER  NX,  NY:  ",$) 

READ ( * , * )  NX, NY 


** 


WRITE ( * , 1100 )  NY0 ,  NY1 
1100  FORMAT (////, 

./ ,  5X,  "********************************************************" 
•/ , 5X, "**  Enter  the  points  which  mark  the  boundaries  of  the  **", 
. / , 5X , " **  parabolic  current  flow.  NY0 ,  NY1  =  ",2I5,3X,  "**", 

./ , 5X, "**  ENTER  NY0 ,  NY1 :  ",$) 

READ ( * , * )  NY0 , NY1 

c  * 

WRITE(*, 1110)  NW,  NT 
1110  FORMAT (////, 

./,5X," ********************************************************  ii 
• / , 5X, "**  Enter  the  number  of  points  for  the  Omega  and  **», 

•/ , 5X, "**  Theta  grid  (spectral  density).  NW,  NT  =  " , 214 , IX"**" , 

. / , 5X , " **  ENTER  NW,  NT:  ",$) 

READ ( * , * )  NW, NT 


** 


WRITE (*, 1120)  NREC 
1120  FORMAT (////, 

./,5 X,"********************************************************n 
•/ , 5X, "**  Enter  the  number  of  spectral  density  files  you  **»| 
. / , 5X , " **  wish  to  create.  NREC  =  "14"  **", 

./, 5X, "**  ENTER  NREC:  ",$) 

READ ( * , * )  NREC 

r  Vfc 

WRITE ( * , 1130)  (NS (N) ,N=1,10) , NREC 
1130  FORMAT (////, 

./, 5X ,  11  ********************************************************  n 
•/ , 5X, "**  Enter  each  Y-location  at  which  spectral  density  **"‘ 
./ , 5X, "**  information  is  desired.  **» 

.  / ,  5X ,  "  **  NS  =  ",10(13)  ,15X,  >.**„' 

• / , 5X , " **  ENTER  all  "14"  values  of  NS:  ",$) 

READ ( * , *)  (NS(N) ,N=1,NREC) 

WRITE ( * , 1140) 
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1140  FORMAT (////, 

./,5X,M ********************************************************", 
./,5X,"**  Enter  the  file  names  for  the  spectral  density  **", 

. / , 5X ,  " **  output  files.  *** 

./,5X," ****************************************************  ****") 
DO  NN  =  1 , NREC 

WRITE ( * , 1150 )  NN,  SPEC(NN) 

1150  FORMAT (/, 

5X ,  M  *  * 

./,5X,"**  File  #  11 , 12 , 11  File  Name  :  f,,A40, 

. / ,  5X,  "**  ENTER  new  filename  -OR-  [/]  : 11  ,$) 

READ ( * , ' (A40) 1 )  NSPEC 

IF  (NSPEC  .NE.  •/')  SPEC (NN)  =  NSPEC 
END  DO 

* 

WRITE (*, 1160) 

1160  FORMAT (////, 

. / , 5X, " ********************************************************11  f 
. / , 5X , n  *  *  Would  you  like  to  change  any  of  the  data?  (y/n)  :  ",$) 

READ ( *  ,  ' (Al) 1 )  ANS 

IF  (ANS  .EQ.  1 Y 1  .OR.  ANS  . EQ .  'y 1 )  GOTO  10 

*  * 


CALL  CLEARED I SPLAY 
WRITE ( * ,  1  ( 5  (/ )  ,5X, 

.  11  ********************************************************  II 

/ 

./,5X,"**  Do  you  wish  to  save  the  data  in  an  input  file  for  **", 
.  /  ,  5X ,  "  *  *  future  use?  (y/n):  ,,,$),) 

READ ( * , 1 (Al) 1 )  ANS 

** 

IF  (ANS  .EQ.  1 Y  *  .OR.  ANS  .EQ.  fy '  )  THEN 
WRITE (*, 1 (5X, 

.  11 ************************************************  ********** 

./,5X,"**  Enter  the  name  of  the  input  file:  M,$)!) 

READ ( * , 1 (A40) f )  NFILE 

OPEN (UNIT—4 0 , FILE=NFILE , F0RM= 1  UNFORMATTED  * ) 

WRITE (40)  TITLE , AKL, AK1 , AK2 , UW , PSIW , UM , U0 , Y0 , Y1 , NX , NY , NY0 , NY1 , 
NW, NT, NREC, (NS (N) ,N=1,NREC) , (SPEC(NN) ,NN=1,NREC) 
CLOSE (UNIT  =  40) 

END  IF 

*  * 

20  UN  =  UM  -  U0 

DY  =  ( Y1  -  Y0 ) /FLOAT (NY  -  1) 

** 

*****  For  a  square  realization  plane. 

*  * 

DO  1=1, NY 

YS ( I )  =  Y0  +  DY* (I  -  1) 

XS(I)  =  DY*(I  -  1) 

END  DO 
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DO  I  =  1, NY 

IF  (I  . LE*  NYO  .OR.  I  .GE.  NY1)  THEN 
U(I)  =  UO 
ELSE 

U(I)  =  ( YS  ( I )  -  YS(NYO)  )*(YS(I)  -  YS  (NY1)  )  *  (-4 . 0*UN/  ( YS  (NY1) 

-  YS (NYO) ) **2)  +  UO 

END  IF 
END  DO 

** 

W1  =  SQRT (GEE*AK1) *UM/GEE 
W2  =  SQRT (GEE*AK2 ) *UM/GEE 
DELTH  =  PI/NT 
YM  =  0.5*(Y1  +  YO) 

WD  =  0.5* ( W2  -  Wl) 

WL  =  SQRT (GEE*AKL) *UM/GEE 

** 

RETURN 

END 

**  ** 

************************************************************************ 
**  ** 

SUBROUTINE  SPECMOD 


**  ** 
************************************************************************ 


*  * 
** 
*  * 
** 
** 
** 
** 
** 
** 
** 
** 
*  * 


This  version  of  SPECMOD  provides  a  realization  for  a  varying 
current. 

NOTE:  XL* *2  has  been  removed  from  the  line  M  821  F(II,J)  =  ...  1 

the  old  instruction  has  been  commented  out.  This  needs 
to  be  done  to  the  other  &INTF@. 


** 
** 
** 
*  * 
** 
** 
** 
** 
** 


This  version  includes  a  Pierson-Moskowitz  spectral  density  and 
has  taken  into  consideration  both  the  incoming  and  outgoing 
waves  by  introducing  additional  terms  to  2ETA  and  an  additional  ** 
random  phase  PH2 . 


** 

** 


************************************************************************* 


**  ** 
COMMON  /INP/  GEE , UM, YO , Y1 , Wl , PI , YM, WD, DY , WL, XS ( 100 )  ,YS(100)  , 

U(100) , NX, NY , PSIW , UW, NW, NT , DELTH , AKL, NREC , NS ( 10 ) , 
TITLE, SPEC (10) 

COMMON  /OUT1/  DELW (99 ) ,THET(99) , THETA ( 99 ) ,WS(99) ,THIR(99) ,THL(99) 
COMMON  /OUT2/  ZETA ( 101 , 101) , ZETA2 ( 101 , 101) , XX ( 101) , YY ( 101) 

COMMON  /OUT3/  F ( 99 , 99 ) , ETA (99 , 99 ) , AK (99 , 99 ) , PH ( 99 , 99 ) , F2 ( 99 , 99 ) 

** 

CHARACTER*40  SPEC 
CHARACTER* 8  0  TITLE 

*  * 

CALL  CLEAR  DISPLAY 
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*  * 

***** 

** 

** 


** 


OPEN (UNIT=66 , FILE= 1  refract . out  ' ) 

Here  we  want  to  expand  the  Omega  and  Theta  grid. 


DO  K  =  1,NW 

DELW(K)  =  2.0*WD/NW 
WS(K)  =  W1  +  0 . 5*DELW (K) * (2*K 
END  DO 


1) 


** 


*  * 

*  *  *  *  * 
***** 
***** 
*  * 


DO  II  =  1,NT+1 

THET(II)  =  -0.5*PI  +  (II  -  1) *DELTH 
END  DO 

DO  II  =  1,NT 

THETA (II)  =  0. 5*  (THET(II-fl)  +  THET(II)) 
END  DO 


This  section  calculates  the  Fourier  coefficients  of  the 
frequency  direction  spectral  density  and  the  coefficients 
of  the  surface  elevation  field  .  F(W(J) ,THET (I) ) . 

TFACTOR  =  GEE/UM 
WRITE ( * , 1 (/////) ') 

DO  1=1, NY 

WRITE  (*,  1  ("»»>  Computing.  Please  wait.  I  =",13)')  I 
XL  =  U(I)/UM 
XLT  =  1.0/XL 
DO  J  =  1 ,  NW 
W  =  WS  (J) 

EEE  =  0 . 5* ( 1 . 0  +  2 . 0*W  -  SQRT(1.0  +  4 . 0*W) ) / (W*W) 

IF  (W  .LT.  0.1)  EEE  =  1.0  -  2.0*W  +  5.0*W*W  -  14.0*W*W*W  + 

42 . 0*W**4  ~  13  2.0  *W*  *  5 

TIR  =  ASIN (EEE) 

AB  =  -XLT/W 

IF  (W  . LT:  XLT)  AB  =  -1.0 
TB  =  -0.5*PI 

IF  (W  .GT.  XLT)  TB  -  ASIN (AB) 

TM  =  0.5* (TIR  +  TB) 

TD  =  0.5* (TIR  -  TB) 

UNIT  «  1.0 

IF  (W  .LT.  WL)  UNIT  =  -1.0 

TL  =  UNIT*  0 . 5  *PI 

STHL  =  XLT* (W  -  WL)/W**2 

IF  (ABS (STHL)  .LE.  1.0)  TL  =  ASIN (STHL) 

THL(K)  =  TL 
DO  II  =  1, NT 

ETA ( II , J)  =  WS(J) *SIN(THETA(II)  ) 

BR  =  4 . 0*XL*ETA(II, J) 

THEX  =  ASIN (-1.0/ (1.0  +  XL*WS ( J ) ) * *2 ) 
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823 

IF  (SIN(THEX)  .GT.  ( -1 . 0/ ( 4 . 0*XL*WS ( J) ) ) )  GOTO  823 

IF  (BR  . LT .  -1.0)  GOTO  820 

GOTO  821 

IF  (THETA (II)  .LT.  THEX)  GOTO  820 

GOTO  821 

820 

F ( II , J)  =  0.0 

GOTO  822 

821 

ETAIIJ  =  2.0*XL*ETA(II,J) 

SQETA  =  SQRT (1.0  +  4 . 0 *XL*ETA ( II , J)  ) 

FIIJ  =  ((SQETA  -  1.0) /ETAIIJ) **4/SQETA*DELW ( J) *WS ( J) *DELTH 

** 

***** 

***** 

***** 

***** 

** 

These  next  instructions  will  insert  the  Pierson-Moskowitz 
incoming  frequency  function. 

S(W)/W  =  0 . 0081 (UM) **6/ (G**4*W**6) *EXP (—0.74 (UM/ (UW*W) )  **4 ) 
where  UW=wind  velocity  in  cm/sec  and  must  be  defined. 

** 

***** 

***** 

***** 

***** 

***** 

** 

F (II, J)  =  FIIJ*0. 0081*UM**6/ ( (WS (J) **6) * (GEE* *4) ) * 

EXP (-0.74* (UM/ (UW*WS ( J) ) ) **4) 

Insert  G (TH)  =  (COS (NU (THINF-PSIW)  )  **2 

where  NU  is  a  parameter  of  the  P/M  model  and  PSIW  is  the  angle 
between  the  wind  and  the  wake  centerline  and  THINF  is  the 
direction  in  the  initial  plane  for  WJ  and  THI . 

We  need  to  insure  that  ABS (NU*THW)  <  PI/2.0 

c 

IF  (F(II,J)  .NE.  0.0)  THEN 

ENU  =  1.0 

ETAIIJ  =  2 . 0*XL*ETA ( II , J) 

SQETA  =  SQRT (1.0  +  4 . 0*XL*ETA (II ,  J)  ) 

SNTHO  =  (1.0  +  ETAIIJ  -  SQETA) / (XL*ETAIIJ*WS (J) ) 

IF  (ABS (SNTHO)  . GE .  1.0)  THEN 

SNTHO  =  (1.0  -  ETAIIJ  + 

5.0* (XL* ETA ( II , J) ) **2) *SIN (THETA (II) ) 

END  IF 

THINF  =  ASIN (SNTHO) 

THW  =  THINF  -  PSIW 

ENUTH  =  THINF  -  PSIW 

c 

THW2  =  PI  -  THINF  +  PSIW 

ENUTH2  =  PI  -  THINF  4-  PSIW 

IF  (ABS (ENUTH2 )  . GE .  0.5*PI)  THEN 

F2 (II , J)  =  0.0 

ELSE 

F2 ( II , J)  =  F (II , J) * (COS (ENUTH2 ) ) **2 

END  IF 

IF  (ABS (ENUTH)  . GE .  0.5*PI)  THEN 

F (II , J)  =  0.0 

ELSE 

F (II , J )  =  F ( II , J) * (COS (ENUTH) ) **2 

END  IF 

END  IF 
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** 

822 

k  k 

kkkkk 
k  k  kk  k 
kk 


kk 


k  k 

kk 

kkkkk 
k  k 


kk 

kkkkk 

kk 


F(II,J)  =  SQRT (2.0*F(II,J) ) *TFACTOR 
F2 (II , J)  =  SQRT ( 2 . 0*F2 ( II , J) ) *TFACTOR 

The  nondimensional  wavenumber,  KIJ  is  now  computed  AK(I,J) 
the  dimensional  K  is  produced  by  multiplying  by  UM**2/G 

ETAIJ  =  4 . 0*XL*ETA(II, J) 

IF  (ETAIJ  .LT.  -1.0)  THEN 
AK ( 1 1 ,  J )  =  0.0 
F ( II , J)  =  0.0 
F2 (II, J)  =0.0 
ELSE 

ETAIIJ  =  2 . 0*XL*ETA (II , J) 

SQETA  =  SQRT (1.0  +  4 . 0*XL*ETA ( II , J) ) 

AK ( II , J)  =  (WS ( J) * (SQETA  -  1.0) /ETAIIJ) * *2 
END  IF 

AKIJ  =  AK ( 1 1 , J ) 

IF  (AKIJ  .LT.  AKL*UM/TFACTOR)  THEN 
AK ( I I , J )  =  0.0 
F ( II , J)  =  0.0 
F2 (II , J)  =  0.0 
END  IF 
END  DO 

IF  (I  .EQ.  NS ( 1 ) )  THEN 

WRITE (66, 1 ( 4X , 12E10 . 5 ) 1 )  (F (K, J) ,K=1, 12) 

WRITE (66, ' ( 4X , 12E10 . 5 ) '  )  (F(K,J) ,K=13,24) 

END  IF 
END  DO 

IF  (I  .EQ.  NS ( 1 ) )  THEN 
DO  J  =  1,25 

WRITE ( 66 , ■ ( 4X , 12E10 . 5 )  f)  (AK(K,J) ,K=1,12) 

WRITE (66, 1 (4X, 12E10.5) 1 )  (AK (K, J) ,K=13,24) 

END  DO 
END  IF 

CALL  SPECTRA_FILE ( I ) 

Initialize  ZETA(II,J) . 

II  =  I 
DO  J  =  1,NY 

ZETA ( J , II )  =  0.0 
ZETA2 ( J , II )  =  0.0 
END  DO 

Random  phase. 

JT  =  1 
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** 

***** 
*  * 


** 

*  *  *  *  * 
** 


** 

*  *  *  ** 
*  **** 
** 


DO  IR  =  1 , NT 
DO  JR  =  1,NW 

PH ( IR, JR)  =  URAND ( JT) *2 . 0*PI 
END  DO 
END  DO 

Define  physical  position. 

COEFF  =  100. 0*TFACTOR/UM 
XX (I)  =  COEFF*XS ( I ) 

DO  JJ  =  1 ,  NY 

YY(JJ)  =  COEFF*YS ( JJ) 

DO  JT  =  1,NW 

Reflection  boundary  THETA  =  THIR. 

WS JT  =  SQRT (1.0  +  4 . 0*WS ( JT) )  -  1.0 

THIR ( JT)  -  ASIN((WSJT/(2.0*WS(JT)*(1.0  -  XL)  +  XL*WSJT) )  **2 ) 
DO  IT  =  1,NT 

These  next  lines  are  needed  to  control  angles  to  be  less  then 
1000*PI  as  requred. 

THP  =  AK ( IT , JT) * (XX (II) *SIN (THETA ( IT) )  + 

YY(JJ) *C0S (THETA (IT) ) )  +  PH(IT,JT) 

THM  =  AK ( IT , JT) * (XX (II) *SIN (THETA (IT) )  - 
YY ( JJ) *COS (THETA ( IT) ) )  +  PH(IT/JT) 

IF  ( ABS (THP)  .GE.  1000. 0*PI)  THEN 
III  =  THP/ ( 2 . 0*PI ) 

THP  =  THP  -  III*2.0*PI 
END  IF 

IF  (ABS (THM)  .GE.  1000. 0*PI)  THEN 
III  =  THM/ ( 2 . 0*PI ) 

,THM  -  THM  -  III*(2.0*PI) 

END  IF 
Y  =  YS(JJ) 

THETAIT  =  THETA (IT) 

THIRJT  =  THIR ( JT) 

IF  (THETAIT  .LT.  THIRJT  .OR.  Y  . GE .  YM)  THEN 
ZETA ( JJ , II )  =  ZETA ( JJ , II )  +  F ( IT , JT) *COS (THP) 

ELSE 

ZETA ( JJf II)  =  ZETA (JJ, II)  + 

F ( IT , JT) * ( COS (THP)  +  COS (THM)) 

END  IF 

THETAIT  =  THETA (IT) 

THIRJT  =  THIR ( JT) 

IF  (THETAIT  .LT.  THIRJT)  THEN 
ZETA2 ( JJ , II )  =  ZETA2 ( JJ , II )  + 

F(IT, JT) * (COS (THP) +C0S (THM) ) 

ELSE 
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** 


** 


ZETA2 ( JJ , II )  =  ZETA2(JJ,II)  + 

F ( IT , JT) *COS (THP)  +  F2 (IT, JT) *COS (THM) 

END  IF 
END  DO 
END  DO 
END  DO 
END  DO 

CALL  PLOT_FILE 

CLOSE (UNIT  =  66) 

RETURN 
END 

** 

********************************************************** *************** 
**  ** 

SUBROUTINE  OUTPUT 

**  ** 

************************************************************************ 
**  ** 

**  This  subroutine  outputs  the  data  files.  ** 

**  ** 

kkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkk 
**  ** 

COMMON  /INP/  GEE , UM , YO , Y1 , W1 , PI , YM, WD , DY , WL, XS ( 100 ) ,YS(100) , 

U(100)  ,NX,NY,PSIW,UW,NW,NT,DELiTH,AKL,NREC,NS(10)  , 
TITLE  f SPEC ( 10 ) 

COMMON  /OUT1/  DELW ( 99 ) ,THET(99) , THETA (99) ,WS(99) ,THIR(99) ,THL(99) 
COMMON  /OUT2/  ZETA ( 101 , 101)  , ZETA2 ( 101 , 101)  , XX ( 101 )  , YY ( 101) 

COMMON  /OUT3 /  F (99 , 99 ) , ETA ( 99 , 99 ) , AK ( 99 , 99 ) , PH (99 , 99 ) , F2 ( 99 , 99 ) 


** 

CHARACTER* 4 0  SPEC 
CHARACTER* 80  TITLE 

** 

k  k  k  -k  k 

Write  the  output 

kk 

ENTRY  SPECTRA  FILE 
NDIM  =  1 

k  k 

DO  NN  -  1,NREC 

*  -k 


IF  (I  .EQ.  NS(NN))  THEN 

OPEN (UNIT=55 , FILE=SPEC (NN) , FORM= ' UNFORMATTED ' ) 

WRITE (55)  TITLE 

WRITE (55)  NT, NW, NDIM 

WRITE ( 55 )  (THETA(II) ,11=1, NT) 

WRITE ( 55 )  (WS (J) , J=1,NW) 

DO  II  =  1 , NT 

WRITE ( 55 )  (F(II, J) , J=1,NW) 

END  DO 
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** 


** 


** 


** 


** 


** 


CLOSE (UNIT  =  55) 
END  IF 
END  DO 
RETURN 

ENTRY  PLOT_FILE 
NDIM  =  1 


OPEN (UNIT 
WRITE (50) 
WRITE (50) 
WRITE (50) 
WRITE (50) 


=50, FILE= ' surf . elv ' , F0RM=' UNFORMATTED ' 
TITLE 

NX, NY, NDIM 
(XX(J) , J=1,NX) 

(YY (II) ,11=1, NY) 


) 


DO  J  =  1 , NX 

WRITE (50)  (ZETA(J,II) ,11=1, NY) 

END  DO 

CLOSE (UNIT  =  50) 

OPEN (UNIT=50 , FILE= ' surf 2 . elv ' , FORM= 1  UNFORMATTED ' ) 
WRITE (50)  TITLE 
WRITE (50)  NX, NY, NDIM 
WRITE (50)  (XX(J) ,J=1,NX) 

WRITE (50)  (YY(II) ,11=1, NY) 

DO  J  =  1 , NX 

WRITE (50)  (ZETA2(J,II) ,11=1, NY) 

END  DO 

CLOSE (UNIT  =  50) 


RETURN 

END 


** 


"k  k 

************************************************** ********************** 
** 

REAL  FUNCTION  URAND(IY) 

** 

*********************************************************^^  ********** 
kk 

k  k 

Urand  is  a  uniform  random  number  generator  based  on  theory  and  ** 

suggestions  given  in  D.E.  Knuth  (1968)  vol.2.  The  integer  IY  ** 

**  should  be  initialized  to  an  arbitrary  integer  prior  to  the  first  ** 

A  call  to  URAND.  The  calling  program  should  not  alter  the  value  ** 

n'P  TV  hofuocm  ^  *1  *1  _  J TTn  k  ATr>  „  _ 


** 

** 

** 


-  - -J  UUL  aitci  1,11 1;  Vdiue 

of  IY  between  successive  calls  to  URAND.  Values  of  URAND  will 
be  returned  in  the  interval  (0,1). 


** 
** 

******************************************************^^*  ^^^^** 
k  k 

SAVE  IA,  IC,  HALFM,  M2,  MIC,  S 

** 
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DOUBLE  PRECISION  HALFM 

** 

DATA  M2 /0/ , ITWO/ 2/ 

** 

*****  If  first  entry,  compute  machine  integer  worlength. 

** 

IF  (M2  .EQ.  0)  THEN 
M  =  1 

10  M2  =  M 

M  =  ITWO*M2 

IF  (M  .GT.  M2)  GOTO  10 
HALFM  =  M2 

*  * 

*****  Compute  multiplier  and  increment  for  linear  congruential  method. 
** 

IA  =  8*IDINT ( HALFM* DATAN ( 1 . 0D0 ) /8.0D0)  +  5 

IC  =  2*IDINT (HALFM* (0.5D0  -  DSQRT ( 3 . 0D0 ) /6 . 0D0 ) )  +  1 

MIC  =  (M2  -  IC)  +  M2 

** 

*****  s  is  the  scale  factor  for  converting  to  floating  point. 

** 

S  =  0 . 5/HALFM 
END  IF 

** 

*****  Compute  the  next  random  number. 

** 

IY  =  IY*IA 

** 

*****  The  following  statement  is  for  computers  which  do  not  allow 
*****  integer  overflow  on  addition. 

** 

IF  (IY  .GT.  MIC)  IY  =  (IY  -  M2)  -  M2 

** 

IY  *  IY  +  IC 

** 

*****  The  following  statement  is  for  computers  where  the 

*****  wordlength  for  addition  is  greater  than  for  multiplication. 

** 

IF  (IY/2  .GT.  M2)  IY  =  (IY  -  M2)  -  M2 

** 

*****  The  following  statement  is  for  computers  where  integer 
*****  overflow  affects  the  sign  bit. 

** 

IF  (IY  .LT.  0)  IY  =  (IY  +  M2)  +  M2 
URAND  =  FLOAT ( IY) *S 

** 

RETURN 

END 

*  * 

********************************************************************** 
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SUBROUTINE  CLEAR  DISPLAY 


** 

*********************************************************.;f*^.fc.^^ 

** 

*****  Home  cursor  and  clear  display. 

** 


*  * 

**** 


CHARACTER* 1  ESC 
CHARACTER* 4  BUFFER 

ESC  =  CHAR (27) 

BUFFER  =  ESC//,H,//ESC//' J' 
WRITE ( * , 1 ( A) ' )  BUFFER (1:4) 
RETURN 
END 


** 

****************************************************** 

** 


************** 


* 

* 


* 

* 


SUBROUTINE  PAUSE_SEC ( ISEC) 

** 

***************************************** 

** 

*****  Pauses  the  program  by  ISEC  seconds 
** 


** 

************  ***************** 


** 


IF  (ISEC 


-LE.  0)  RETURN 


WRITE  (* ,  1  (8  (/)  ,  ">  Type  <  /  >  to  continue  when  ready:  "  ,$)f) 
READ ( * , * )  SLASH 

RETURN 

END 
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